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FOREWORD 


This report. was prepared by Universal Analytics, Inc. (UAI) under NASA 
Contract NASI -11677, for the National Aeronautics and Space Administra- 
tion, Langley Research Center, Hampton, Virginia. Dr. E.I. Field of UAI 
served as Program Manager with Mr. S.E. Johnson contributing significantly 
to the program's success. This work was administered under the direction 
of NASTRAN Systems Management Of f ice . (NSMO) with Dr. J.P. Raney and Mr. H. 
Adelman as the contract and technical monitors, respectively. UAI wishes 
to acknowledge their continued support. and encouragement provided through- 
out the project. 

This report is presented in one volume, containing unclassified material. 

The objective of this contract effort. was to incorporate the three-dimen- 
sional family (linear, quadratic . and . cubic) of isoparametric solid elements 
into NASTRAN, NAsa's STRu ctural ANalysis computer program. The program 
code was designed. and initially written for eventual installation into 
Level 16 of NASTRAN which is still in development. Actual implementation, 
however, was carried out in the Level 15.1 version that was made available 
by NSMO on the UNIVAC 1108 at the Jet Propulsion Laboratories in Pasadena, 
California. 

The completed NASTRAN modifications and additions were successfully demon- 
strated by exercising three of the Rigid Formats for static, normal modes 
and buckling analysis. Final program documentation was prepared in the 
form- of new sections and updates for direct insertion into the NASTRAN 
Theoretical, User's and Programmer's Manuals. These inserts are incorpor- 
ated in the Appendices of this report. Compilation listings, program decks 
and the computer printout of the demonstration problems, submitted under 
separate cover, complete the deliverables for this contract. 
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ABSTRACT 


This report summarizes the work done by Universal Analytics, Inc. for 
NASA, Langley Research Center, to implement the three-dimensional family 
of linear, quadratic and cubic isoparametric solid elements into the 
NASA Structural Analysis program, NASTRAN. This work included program 
development, installation, testing and . documentation. The addition of 
these elements to NASTRAN provides a significant increase in modeling 
capability particularly for structures requiring specification of temper- 
atures, material properties, displacements and stresses which vary through- 
out each individual element. 

Complete program documentation is presented in the form of new sections 
and updates for direct insertion to the three NASTRAN manuals. The results 
of demonstration test problems are summarized. Excellent results are 
obtained with the isoparametric elements for static, normal mode and buck- 
ling analyses . 
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INTRODUCTION 


The NAsa STRuctural ANalysis program (NASTRAN)-. has been developed with the . 
intent of being responsive to a wide variety of user needs in the area of 
matrix structural analysis. In order to ensure the future usefulness of 
NASTRAN, it is essential that the best state-of-the-art finite element 
technology be represented in the system. New aerospace vehicle concepts, 
such as the Space Shuttle, have given impetus to NASA, through Langley 
Research Center (LRC) , for undertaking the task of updating and improving 
NASTRAN. In response to this need, Universal. Analytics , Inc. (UAI) has 
developed and installed the three-dimensional family of linear, quadratic 
and cubic isoparametric solid elements into NASTRAN. With this effort now 
completed, these three new elements greatly enhance NASTRAN 's capability to 
solve any three-dimensional solid problem such as would be involved with both 
the structural and thermal design and analysis of the Space Shuttle. 

The complete documentation . of the new isoparametric solid element capability 
has been prepared in the form of new sections and updates to be inserted 
directly into each of the NASTRAN Theoretical, User's and Programmer's 
Manuals. These inserts are presented in this report in Appendixes A, B, 
and C respectively. The next section of this report discusses these docu- 
mentation inserts as they were derived from the original program specification. 

The remainder of this report describes the program implementation, testing 
and demonstration. Though the program was designed and coded for eventual 
implementation into Level 16 of NASTRAN when it becomes available, this 
report also documents the actual implementation as carried out in Level 15.1. 

A tabulation of typical computer run times is included for purposes of 
comparison. Instructions to the NASTRAN Maintenance Contractor for implemen- 
tations of the new routines into Level 16 are presented in Appendix D. The 
excellent results obtained from the demonstration test problems are summarized 
from the reports presented in Appendix E. The observations made from testing 
the new element capability are presented in the concluding section of this 
report. 
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DEVELOPMENT OF SPECIFICATIONS 


The complete final documentation for the installation of the three isopara- 
metric elements into NASTRAN is presented in the appendices of this report. 
The formal documentation of the theory, user input and program additions 
and modifications was prepared for direct insertion into the corresponding 
three NASTRAN Manuals. These inserts include new sections and appropriate 
up dates to the existing sections. The set of inserts for each manual 
is described in the paragraphs below. 

The original specifications, consisting of both Theoretical and Functional 
Module Mathematical Specifications (FMMS) , were prepared for installation 
of the new elements into Level 16 of NASTRAN. Level 16, when completed, 
will provide direct access capabilities required for efficient assembly of 
the isoparametric element matrices. Actual installation and testing, however 
had to be carried out in Level 15.1, the latest available version of NASTRAN. 
Instructions to. the NASTRAN maintenance contractor for installation into 
Level 16 when available are presented in Appendix D. 

Theoretical . Specifications 

The documentation of the theory for isoparametric elements is presented in 
Appendix A along with instructions for their insertion into the NASTRAN 
Theoretical Manual . . A complete mathematical description is provided with 
references. Both the structural. and heat transfer theory is presented 
showing the distinct advantages gained by the direct integration approach 
to generating the element matrices. The unique advantage of the isopara- 
metric family of elements is shown to derive from the assumption that temper- 
atures, material properties, displacements . and stresses may vary throughout 
the element. That is, these elements more accurately represent the structure 
and its response than do conventional constant strain elements. 

User's Specifications 

The inserts and updates for the NASTRAN. User ' s Manual are presented in 
Appendix B. These include the descriptions of new Bulk Data cards 
(CIHEX1, CIHEX2 , CIHEX3 , PIHEX and PL0AD3) , changes to the description of 
old Bulk Data cards (TEMP and TEMPD) , additions to the acceptable element 
plot "SET" and associated symbols, and instructions to the user for modeling 
with this new capability. The new error messages, which may be issued during 
processing, are shown as additions to the User's Manual error message list. 

Program. Specifications 

By far the bulk of the detail documentation involves updates and additions 
to the NASTRAN Programmer's Manual . These manual inserts are incorporated 
in Appendix C. They include all modifications to existing NASTRAN routines 
and block data descriptions as well as the description of all new subroutines 
The specifications for the matrix operations and sequencing of integration 
steps required to implement the theoretical specifications are included as 
updates to the sections describing the appropriate NASTRAN modules which 
perform these tasks. Instructions to the NASTRAN maintenance contractor 
for installation into Level 16 are presented in Appendix D. The organization 
and development of the Level 16 routines. and their adaptation for testing 
in Level 15.1 is summarized in the following section. 
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PROGRAM IMPLEMENTATION 


This section describes the development and implementation of the family of 
isoparametric solid elements for NASTRAN. The design of subroutines for 
Level 16 which was then modified to function for Level 15 is discussed. 
Modifications which were made to functional modules are described. The 
design and development of new subroutines are also described. Historical 
execution times are presented for an approximate comparison on computa- 
tional efficiency for all three elements. The complete documentation of 
the program changes and additions are presented in Appendix C as inserts 
and updates to the NASTRAN Programmer's Manual. 


Program Development 

The element matrix generation subroutine XIHEX was initially developed in 
a stand-alone driver. It was designed to be installed Level 16. Sub- 
routines KIHEX, MIHEX, and DIHEX were extracted from XIHEX to generate the 
element matrices in the Level 15 environment. These subroutines were 
developed as an interim measure so that the new elements could operate in 
Level 15. As KIHEX, MIHEX, and DIHEX were debugged, all changes and 
corrections made in them were also made in XIHEX. However, XIHEX could not 
be completely checked out because it could not operate in Level 15 and the 
stand-alone driver was inadequate for the task. Thus, it must be checked 
when it is installed in Level 16. Instructions and suggestions for this 
installation are contained in Appendix D. They should be helpful to the 
maintenance contractor in pointing out possible problem areas and incomplete 
interfaces with new Level 16 subroutines. Subroutine QIHEX, which computes 
the heat generation load vector for the isoparametric solid elements, could 
not be installed in Level 15. It was designed and developed for Level 16. 
Appendix D also describes the implementation of this subroutine. All other 
subroutines modified or added. were developed. in the environment of NASTRAN 
Level 15 and no significant changes are necessary to install them in Level 16. 


Program . Modifications 

Many existing NASTRAN functional modules and subroutines were modified to 
implement the isoparametric solid elements. Several new subroutines were 
also added. These modules and subroutines are listed in Table 1 and brief 
descriptions of the changes to each functional module are described in 
Table 2. Many of the changes are those normally required when implementing 
new elements. However, in this case, changes were also required in the PL0T 
module (to plot 3-dimensional elements) , the GP3 module (to process a new 
external load) , and in various other modules to accommodate the large space 
requirements of the cubic element. All modifications to existing subroutines 
and new subroutines contain extensive in-program documentation with comment 
cards to assist in any future program development. Further documentation is 
located in the inserts and updates for the Programmer's Manual, Appendix C. 

These inserts and updates are based on the Programmer's Manual for Level 14 
since the Level. 15 Programmer !s Manual. was not completed . during the fulfillment} 
of this contract. 
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TABLE 1. SUBROUTINES MODIFIED OR ADDED TO IMPLEMENT THE FAMILY 
OF ISOPARAMETRIC SOLID ELEMENTS 


Module 

Subroutine 

Module 

Subroutine 

XCSA 

LD01 

GP2 

GP2 


LD02 




LD03 

PLTSET 

CNSTRC 


LD04 


C0MECT 


LD05 

PL0T 

ELELBL 


LD06 


SHAPE 


LD07 




LD08 

GP3 

GP3 


LD09 


GP3A 


LD10 


GP3BD 


LD11 


GP3C 


LD12 


GP3D 


LD13 

TA1 

TA1A 

IFP 

IFP 


TA1B 


IFS5P 


TA1ETD* 


IFX1BD 




IFX2BD 

SMA1 

KIHEX* 


IFX3BD 


SMA1A 


IFX4BD 


SMA1BD 


IFX5BD 

SMA2 

MIHEX* 


IFX6BD 


SMA2 


IFX7BD 


SMA2A 




SMA2BD 


Module Subroutine Module Subroutine 

SSG1 EDTL Various GPTABD 
EXTERN IHEXSD* 

IHEX* IHEXSS* 

PERMUT> PREMAT 

PL0ADS* TKTZTK* 

QIHEX** 

EMG XIHEX** 

DSMG1 DIHEX* 

DS1 

DS1A 

DS1ETD* 

SDR2 SDRETD 

SDR2B 
SDR2BD 
SDR2D 
SDR2E 
SIHEXl* 

SIHEX2* 

0FP 0FP 

0FP1A 
0FP1BD 
0FP5BD 
0F1PBD 
0F2PBD 
0F3PBD 
0F4PBD 


New subroutine. 

**New subroutine not implemented in Level 15. 
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TABLE 2. DESCRIPTIONS OF FUNCTIONAL MODULE MODIFICATIONS 


XCSA - New entries were made in the card name table in sub- 
routines LDi’, i=l, 13 for the new: Bulk Data cards 

IFP - New Bulk Data cards were added 

GP2 - Array sizes were increased to accommodate elements 
with 32 grid points 

PLTSET - Array sizes were increased to accommodate elements 
with 32 grid points 

PL0T - Capability for plotting solid elements was implemented 

GP3 - PL0AD3 card processing was implemented 

TA1 - Capability to append grid point temperatures to EST/ECPT 
entries was implemented 

SMA1 - Stiffness and conductance matrix generation for the new 
elements was implemented 

SMA2 - Mass and capacitance matrix generation for the new 
elements was implemented 

SSG1 - Load vector generation for thermal and pressure loads 
on the new elements was implemented 

DSMG1 - Differential stiffness matrix generation for the new 
elements was implemented 

SDR2 - Stress calculations for the new elements was implemented 


0FP 


Stress printout formats for the new elements were 
implemented 




Program Additions 


. The New NASTRAN subroutines required to implement the family of isoparametric 
J solid elements have been designed. to operate in the Level 16 environment. 
vThe subroutine (XIHEX) which generates stiffness, mass, differential stiffness, 
conductance, and capacitance matrices was specifically designed for Level 16 
Functional Module EMG (Element Matrix Generator) . It was altered to create 
four interim subroutines (KIHEX, MIHEX, DIHEX, TKTZTK)/ which will operate in 
Level 15. When the isoparametric solids are incorporated into Level 16, 
these four subroutines should be discarded . and subroutine XIHEX would take 
their place. 

Subroutine XIHEX is divided into three major sections. These are initialization 
and geometry verification, matrix generation, and transformation and output. 

The first section initializes. variables for the matrix generation process and 
checks the geometry of the element. If the geometry is acceptable, matrix 
generation proceeds. UAI believes that a good program design would verify 
the geometry of all . elements , particularly isoparametric elements, before 
generating any element matrices. This is not possible, however in the current 
Level 15 design. An error in the geometry of the last element would result 
in the loss of all element matrices previously computed. For problems involv- 
ing many isoparametric elements, this loss could be quite expensive. 

In the second section of XIHEX, the stiffness (or conductance) and mass (or 
capacitance) and the differential stiffness matrices are generated using •;> 
numerical integration. Only the lower triangular portion of these matrices 
is generated and the matrix is generated. in the basic coordinate system. The 
stiffness and mass (or conductance and capacitance) matrices may be generated 
in the same call to subroutine XIHEX, or they, may be generated in individual 
calls. If they are being generated in the same call, they are generated simul- 
taneously or serially depending on the amount of open core available. 

In the third section of XIHEX, the symmetric halves of the mass and stiffness 
matrices are first transformed . to the global, coordinate system. Then-.three 
complete . columns , corresponding to one grid. point, are extracted from the 
matrices and passed. to subroutine EMG0UT for output. Conductance and capaci- 
tance matrices are transposed in core and the full matrix is passed to EMG0UT. 

Subroutines . KIHEX, MIHEX, and DIHEX are interim subroutines for use in Level 15. 
They perform the same operations as XIHEX. KIHEX computes stiffness and conduc- 
, tance matrices; MIHEX computes . mass and capacitance matrices; and DIHEX computes 
the differential. stiffness matrix.. These subroutines call the interim sub- 
routine TKTZTK, which was developed. to avoid a compiler error on the UNIVAC 1108. 
This error was avoided in XIHEX whic h, therefore, does not call TKTZTK. KIHEX, 
MIHEX, and DIHEX were extracted fromIXIHEX. To improve computational efficiency, 
they generate only those portions of the element matrix required for each pivot 
point. Also, because materials with thermal capacity are not available in Level 
15, MIHEX uses the value of the thermal conductivity coefficient to emulate the 
element "capacitance" matrix. 
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Two new subroutines were required . to implement the thermal and pressure 
load for the isoparametric solid elements. Subroutine IHEX generates 
the thermal load vector and subroutine PL0AD3 ; generates the pressure load 
vector. Following NASTRAN programming standards, both were originally 
coded in single precision. However, execution of the symmetrical demon- 
stration problems produced noticeably unsymmetric results. This problem 
was traced to the load vector generation, where small extraneous numbers 
were being generated in elements of the load vector which should have been 
zero. IHEX and PL0AD3 were then modified to operate in full double pre- 
cision, truncating only the final element load vector on being added to the 
global load vector. This decreased the size. of the "small extraneous 
numbers" by several orders of magnitude and. greatly improved the symmetry 
of the results. The increased accuracy is attributed mainly to the more 
accurate computation of the isoparametric. shape functions, their derivatives, 
and of the Jacobian matrix. 

Several other new subroutines were also required. IHEXSD and IHEXSS compute 
the shape functions, their derivatives, and the determinant and inverse of the 
Jacobian matrix in double or single precision. SIHEX1 computes stress 
matrices and SIHEX2 computes stresses for the new elements. QIHEX calculates 
the heat generation load vector for. a unit value of external heat per unit 
volume. QIHEX is not implemented in Level 15. Appendix D describes how it 
may be implemented in Level 16. 

Program Execution 

The execution times required to generate . element matrices in Level 16 
using subroutine XIHEX will be substantially less than the execution times 
of Level 15 . Table 3 shows the approximate matrix generation times for 
subroutine XIHEX operating with the stand-alone driver used in its initial 
development. These execution times are considerably less than those of 
Level 15 which are presented in Table 4. This is so even though the Level 
15 times include the overhead of executing a complete functional module and 
not just a single subroutine. 

For complete problem execution. Table 5 presents approximate computer run 
times for each of three Rigid Formats (1 - Static, 3 - Normal Modes, and 
5 - Buckling). These timings were for the three demonstration problems 
described in the following section under Testing and Demonstration. For 
both Rigid Formats 3 and 5, the inverse power method was used.' For purposes 
of comparison, the characteristic dimensions of the problem are also shown 
in Table 5. 

Of noteworthy interest is the comparison of run times for the relative 
accuracy obtained in the solutions. That is, the IHEX3 model in every case 
has the shortest run time and yielded the most accurate results when deflec- 
tions were used for comparison. Comparing element stresses, however, the 
IHEX2 model yielded the best results for the thermal loading conditions. 

The next section of this report presents the details of the problems 
analyzed. 
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TABLE 3. EXECUTION TIME FOR THE LEVEL 16 ELEMENT MATRIX- 
GENERATION SUBROUTINE XIHEX 


Matrix 

Type 

Element 

Type 

Number of 
Integration 
Points 

UNI VAC 1108 
CPU Time 
(Minutes) 

Stiffness 

IHEX1 

2 

0.005 


IHEX2 

3 

0.1 


IHEX3 

4 

0.6 

Mass 

IHEX1 

2 

0.002 


IHEX2 

3 

0.05 


IHEX3 

4 

0.3 


TABLE 4. MODULE EXECUTION TIMES FOR A ONE-ELEMENT MODEL 

IN NASTRAN LEVEL 15 


Matrix 

Type 

Module 

Name 

Element 

Type 

Number of. 
Integration 
Points 

UNIVAC 1108 
CPU Time 
(Minutes) 

Stiffness 

SMA1 

IHEX1 

2 

0.02 



IHEX2 

3 

0.3 



IHEX3 

3 

0.7 



IHEX3 

4 

2.1 

Mass 

SMA2 

IHEX1 

2 

0.01 



IHEX2 

3 

0.1 



IHEX3 

3 

0.3 

Differential 

DSMG1 

IHEX1 

2 

0.05 

Stiffness 


IHEX2 

2 

0.2 



IHEX3 

2 

0.5 
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TABLE 5. SOLUTION EXECUTION TIMES FOR DEMONSTRATION PROBLEMS 
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TESTING AND DEMONSTRATION 


Reports on the final demonstration testing for static normal mode and 
buckling analyses are incorporated in Appendix E. The demonstration 
problem selected and the results obtained for each model are discussed. 

At each step in the program development, detail numerical checking was 
performed to qualify the logic and the computations. The routines for 
generating the element matrices were then integrated and exercised in 
a skeleton program to simulate the data blocks and control parameters 
to be passed from one NASTRAN module to another. At this level of 
development only the element matrix generation routines were tested 
to verify their design for implementation into Level 16 of NASTRAN. 

These routines were then modified as described in the preceding section 
for insertion into Level 15. The necessary modifications were made to 
existing NASTRAN modules to prepare for their insertion and again each 
subroutine was tested for logic and correct computations. Whenever 
possible, existing programs previously developed by UAI were used to 
verify the final output matrices. Detail hand calculations for both 
structural and heat transfer computations were performed at each inter- 
mediate processing step. Checks were made to verify the shape functions, 
the Jacobian and its inverse, the computations of normals, the transfor- 
mations, the stresses at each integration point, the strain matrix and 
the summation over each integration point. The final element structural 
matrices were then used to solve single element test problems to check 
the overall flow of all other module modifications, loads generation, 
plotting and stress computations. First, simple conditions were specified 
to check equilibrium. Then specific transformations were applied and the 
results again were checked. Combinations of a few elements were then used 
to verify the assembly of element structural matrices and to exercise 
the Rigid Formats to be used for final demonstration problem testing. The 
final demonstration testing was then performed, the results of which are 
summarized below. 


Demonstration Problems 


A cantilever beam was chosen for testing all three isoparametric solid 
elements in NASTRAN. This problem was chosen for two reasons, detail 
theoretical solutions are well known and solid finite elements character- 
istically do not perform well when used to model structures which deform 
primarily in a bending mode. 

Three models were prepared, one with each of the three elements: CIHEX1, 

CIHEX2 and CIHEX3. All three beam models were of uniform crossection 
measuring 12 x 24 x 144 as pictured below. 
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Four static loading conditions were imposed: 

a. Linear thermal gradient through its depth 

(T = 120° @ Y = 0, T = - 120° @ Y = 24) 

b. Uniform temperature rise (AT = 100°) 

c. Compressive axial pressure (used also for buckling) 

(P = 42,837 @ Z = 144) 

d. Uniform . transverse pressure 

(P = 100 @ Y = 0) 

y 

The total number of. elements and grid points used, and the half -band width 
for each model is presented in the following table: 


Model Type 

Number of 
Elements 

Number of 
Grid Points 

Half-Band Width 
(Grid Points) 

Model 1 - IHEX1 Elements 

216 

364 

34 

Model 2 - IHEX2 Elements 

36 

275 

52 

Model 3 - IHEX3 Elements 

8 

148 

44 
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Examination of initial results . which showed lack of symmetry, suggested that 
the load generation routines be changed to double precision, a departure from 
NASTRAN standard conventions. The resulting load vectors, when truncated to 
single precision for storage, provided excellent results. The experience 
gained with the thermal gradient loading condition confirmed that the selection 
of number ^of integration points (specified on the PIHEX card) is critical to 
the accurate computation. of _ thermal loading .vectors. The user . is cautioned to 
request at least three integration points (NIP = 3 on the PIHEX card) for linear 
variations in temperature over the element. No distribution of temperature 
more severe than a "gentle" quadratic should be used. In this case at least 
four integration points are required (NIP = 4) . 

Typical computation times for generation and assembly of the isoparametric 
element matrices are presented in the preceding section. A summary of 
results reported in Appendix E for all three models is presented below for 
each of the three NASTRAN Rigid Formats used. 


Static Analysis (Rigid Format - 1) 

All four static loading conditions described earlier were, imposed on each 
of the three models. The following table presents the results comparing 
the deflection at the tip to the theoretical solution. 


Loading Case 
Maximum Displacement 

NASTRAN 

Theoretical 

Solution 

Model No. 1 
IHEX1 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

Y-Disp. - Load Case a 

1.444 

1.548 

1.533 

1.481 

Z-Disp. - Load Case b 

.2113 

.2104 

,2088 

.2056 

Z-Disp. - Load Case c 

-.2039 

-.2042 

-.2047 

-.2056 

Y-Disp. - Load Case d 

. 1422 

.1561 

.1569 

.1586 


Note that the thermal loading . results . fall on either side of the theoretical 
and the mechanical loading results fall below theoretical, as would be expected. 


Normal Modes (Rigid Format - 3) 

The first three normal modes were extracted for each of the three models using 
the inverse power method and the coupled mass matrix generated by the 
isoparametric element routines. The models were too large to use the Givens 
method, through on smaller test cases, excellent comparisons were obtained. 

The natural frequencies obtained are compared with the theoretical solution in 
the following table. 
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NASTRAN 


Mode 

No. 

Description 

Model No. 1 
IHEX1 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

Theoretical 

Solution 

(cps) 



Freq. 

1 

% 

Error 

Freq. 

% 

Error 

Freq. 

% 

Error 

E 

First Bending Mode 
in the X-Direction 

22.0 

+18 

18.6 

| 

0 

18.6 

0 

18.6 

2 

First Bending Mode 
in the Y-Direction 

38.3 

+4 

36.5 

-2 

36.5 

-2 

37.3 

3 

Second Bending 
Mode in the 
X-Direction 

135.3 

+16 

114.3 

-2 

113.3 

-3 

116.8 


These results show excellent comparisons for the two higher order element 
models. The linear element model yields good results for the second mode 
only. However, because there were three elements through the thickness in 
the x-direction and six elements in the Y-direction, these results are as 
expected. They only point to the well known fact that many linear elements 
are required to adequately model bending behavior. 

Buckling Analysis (Rigid Format - 5) 

The final test for all three models was to compute the critical axial 
pressure buckling load (static loading c) . Three buckling modes were 
computed using the inverse power method. These results are summarized 
below in comparison with the theoretical solution. 


Mode 

No. 

Description 

NASTRAN 

Theoretical 

Solution 

X 

Model No. 1 
IHEX1 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

X 

% 

Error 

m 

% 

Error 

X 

% 

Error 

■ 

X-Direction 

1.406 

40.6 

1.002 

H 

1.001 

.1 

1.0 

2 

Y-Direction 

4.391 

9.8 

3.981 


3.979 

.5 

4.0 

3 

X-Direction 

12.809 

42.3 

9.037 

D 

8.934 

.7 

9.0 
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Again, excellent results are obtained for both the higher order element 
models. And again, the linear element shows much better results in bend- 
ing in the Y-direction as expected. This further emphasizes the need for 
using many linear elements to adequately model bending behavior. 


SUMMARY AND CONCLUSIONS 


Three new solid elements have been implemented into NASTRAN for static, 
dynamic, buckling and heat transfer analysis. These are the linear, quad- 
ratic and cubic members of the isoparametric family of three-dimensional 
elements. The element matrices for both structural and heat transfer 
analyses are generated by direct integration over the element. These 
isoparametric elements have the unique advantage that variations in temper- 
atures, material properties, displacements and stresses may be specified 
throughout the element. These elements can therefore represent a structure 
and it response more precisely than can the conventional constant strain 
elements . 

The basic element matrix generation routines for these elements were designed 
and coded for implementation into NASTRAN Level 16 when it becomes avail- 
able. For testing purposes, however, these routines were modified for and 
implemented into the latest available version of NASTRAN, Level 15:1' 1 Final 
complete documentation, presented in Appendices A-C, was prepared for direct 
insertion into each of the NASTRAN Theoretical User's and Programmer's 
Manuals. Instructions are given in Appendix D. to assist on the eventual . 
implementation into Level 16. Reports are included in Appendix E on testing 
results obtained from demonstration problems exercising three Rigid Formats 
for static, normal modes and buckling analyses. 

Complete checking at each step of the computations was performed manually and 
by comparison with other existing programs. Final verification of the inte- 
grated system was performed using Level 15. For each of the three NASTRAN 
Rigid Formats demonstrated, the results compared very favorable with the theo^ 
retical solution. 

All three elements produced good demonstration test results for static, normal 
mode and buckling analyses. Problems for which thermal gradients were 
specified through the element yielded excellent results when sufficient integra- 
tion points were specified for computation of the thermal load vector. As 
expected, the linear element results showed that it is best used when primarily 
a shearing type behavior is anticipated. And finally, the contention that 
the isoparametric three-dimensional solid elements are superior was validated 
by the excellent test results obtained with the quadratic and cubic elements 
when used to model the primarily bending type behavior of a cantilever beam. 

The implementation of these three isoparametric solid elements, therefore, 
does greatly enhance the total modeling capability of NASTRAN. 
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APPENDIX A 

THEORETICAL MANUAL INSERTS 


The following sections contained in this appendix should be 
added to the NASTRAN Theoretical Manual in their entirety: 

5.13 . ISOPARAMETRIC HEXAHEDRON SOLID ELEMENTS 

7.4 ISOPARAMETRIC HEXAHEDRON SOLID ELEMENTS 

8.2.2 Three-Dimensional Isoparametric Solid 
Element Heat Transfer Matrices 


The remaining sections are marked to indicate changes to the 
existing TM document. 
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5.13 ISOPARAMETRIC HEXAHEDRON SOLID ELEMENTS 

Hexahedron solid isoparametric elements may be used to analyze any three- 
dimensional continuum composed of isotropic or anisotropic materials. 

Examples include thick inserts in rocket engine nozzles, thermal protection 
system insulations, soil structure interaction problems, and geometrically 
complex thick-walled mechanical components such as pumps, valves, etc. The 
isoparametric solid elements take into account: 

1. Isotropic or anisotropic temperature-dependent material properties 

2. Pressure and temperature loads 

3. Coupled mass matrix 

Although solid elements employ only three degrees of freedom at each grid 
point '(the three displacement components), they may be combined with all 
other NASTRAN elements. 

The isoparametric solid elements were first presented by Irons, Ergatoudis 
and Zienkiewicz [Refs. 1, 2, 3 & 4]. They are also called conformal higher 
order elements, since the displacement of the element can be represented to 
any degree one desires, and still maintain interelement compatibility, by 
using more grid points per edge to define element geometry and deformation. 

In practice, however, isoparametric solid elements employing either eight, 
twenty or thirty- two grid points have been found to be adequate to solve most 
problems (Figure 1) . These elements correspond to assuming a linear, para- 
bolic and cubic variation of displacement, respectively. Clough [Ref. 5] 
conducted an evaluation of three-dimensional solid elements and showed that 
the isoparametric elements were superior to other solid elements. He further 
pointed out that the choice of which isoparametric element is best to use 
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\TH 


(a) Linear 



depends on the type of problem being solved. For problems involving plate 
bending type deformations, the higher order elements appear to be best, while 
the linear element is recommended for problems in which shear stresses are 
likely to be large. It is for this reason that all three isoparametric ele- 
ments have been incorporated into NASTRAN. 

The isoparametric elements governing equations are based on minimum energy 
principles. The derivation of these equations consists of assuming a dis- 
placement function for the element which depends on its grid point displace- 
ments, substitute these displacements into the Potential Energy, and minimize 
the energy functional to obtain the governing equations. The detailed 
derivation is next presented. 

5.13.1 Displacement Functions for the Isoparametric Element 
The name isoparametric comes from the fact that the same interpolating func- 
tions are used to represent both the geometry and the deformation of the 
element. This choice insures that the element displacement functions satisfy 
the criteria necessary for convergence of the finite element analysis [Ref. 4]. 
Referring to the curvilinear coordinates (£,r|,£) shown in Figure 1, the 
(x,y , z) coordinates at any point in the element are obtained from the "basic" 
rectangular coordinates at each of the NGP grid points by: 


r 

i- 

I n&p 

'*) 

f 


V 

{ fr. 

~ 1 



where N^C^n,?) are interpolating functions which depend on the number of 
grid points used to define the element geometry. The functions are 
either linear, parabolic or cubic, and correspond to employing two, three or 
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four grid points, respectively, along each edge of the element. This choice 
insures that there are no geometric gaps between grid points. Expressions 
for the interpolating functions and their derivatives are presented in Table 1. 


Now, we represent the deformation of the elements with the identical interpo- 
lating functions used to define the geometry, i.e.: 



NG-P 

E w i(€.n.c) 


-L- 1 



( 2 ) 


where u, v and w are displacements along the x, y and z basic coordinate axes, 
and {u g } represents the vector of grid point displacements. The displacement 
functions, N^(^, r| , £) satisfy the required convergence criteria of adequately 
representing a constant strain state, and insure interelement compatibility 
along the complete element boundary [Ref. 4] 


5 • 13 . 2 Strain-Displacement Relationship 

The strains at any point within the element are given by the well known 



(3) 


Substituting Eq. (2) into (3) yields the relationship for the strain vector 


in terms of the grid point displacements: 




I 
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(4) 
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TABLE la. ISOPARAMETRIC SHAPE FUNCTIONS AND THEIR DERIVATIVES FOR LINEAR ELEMENT - 8 GRID POINTS 
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TABLE lc. ISOPARAMETRIC SHAPE FUNCTIONS AND THEIR DERIVATIVES FOR CUBIC ELEMENT - 32 GRID POINTS 
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In order to evaluate the strain matrr 53 [C], v/e'must calculate the derivatives 
of the shape functions, N^, with respect to x, y and z. Since is defined 
in terms of £, r) and £, it is necessary to use the relation that: 


[J] 


9 to± X 

[ 


where [J] is the Jacobian matrix and is easily evaluated by noting that: 

f irj ail fah/, 2M, F . 1 
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|_ J |_ J (_ ^N&pJ 

The derivatives of the shape functions with respect to p and t; are given 
in Table 1 , and x^, y^ and z-^ are the coordinates of the element grid points, 
and NGP is the number of element grid points. 


5.13.3' Stress-Strain Relations 

The stress-strain relations for a general elastic anisotropic material are: 
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and {a} is the stress vector in the basic rectangular Cartesian coordinate 
system, [G e ] is, in general, a full symmetric material elastic modulus 
matrix, e is the total strain vector 'given by Eq-. (3) , and e is the thermal 
strain: 
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' «*] = {^e} T 


where {a e } is a vector of thermal expansion coefficients, and T is the tem- 
perature distribution within the element, and is determined from user-specified 
grid point temperatures, {T e } by applying the interpolating relation: 

' f = I 5)T.c = [n] K) (10) 

In NASTRAN, the user may specify either isotropic or anisotropic material 

properties. Furthermore, for anisotropic materials, the user may specify the 

properties with respect to a particular orientation that does not necessarily 

coincide with the basic rectangular coordinate system. NASTRAN will then 

transform the material properties to the basic coordinate system by construct- 

M 

ing a transformation matrix, T , such that the material properties in the 
basic coordinate system are given by: 



M • bf K] M ] 

k 

N-MW . 


( 11 ) 
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where [G m ] and {a m ) are the material modulus matrix and thermal expansion 
coefficient vector specified by the user. 


5.13.4 Stiffness, Mass and Load Matrices 

The stiffness, mass and load matrices for the isoparametric element may now 
be derived by application of the Virtual Work Principle: 


6 V - S W = o 


( 12 ) 


where 6U is the internal strain energy in the element due to a virtual dis- 
placement, Su, and 6W is the work performed by the external loads during the 
virtual displacement, i.e.: 

sv = J v H t M dV 

8 W = f [ <5 a.] | d.5 4- \ ^ ^ ^ 

5 J V yj 



where {p} and {f} are 3x1 vector representing surface pressure and body 
forces, respectively, in the x, y and z directions, and the integrations are 
performed over the element volume, V, and surface area, S, on which the pres- 
sure load p acts. 


Substituting Eqs. (2), (4) and (8) into (13) and applying the Virtual Work 
Principle, Eq. (12), yields: 



(14) 


or for any virtual displacement, we have: 



( 15 ) 
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( 16 ) 


where [Kgg ] is the element stiffness matrix and (F e ) is the element load 
vector to surface pressures and temperature. 


Stiffness Matrix 





L & e ] [ C ] A V 


where the infinitesimal volume is in terms of the curvilinear 
coordinates, £, r) and £, 

dV = A^ A* = Jut. £jJ A| A£ 


(17) 


Surface Pressure Load 

r T 


{■>] - J M M A5 


Thermal Load 


T 

j 


V 


t c ] [^e] f^e ' M 


T, 


d.V 


(18) 


(19) 


For dynamic problems, the mass matrix is also required and is easily derived 
by adding the kinetic energy to Eq. (12). The result is: 


Mass Matrix 

f T 

[m] - [n] [nJ pav 

J v 

where p is the mass density. 


( 20 ) 
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;5 - 1 3 . 5 . Numerical Integration 

The integrals in the isoparametric element stiffness, mass and load matrices 
are evaluated by the use of numerical integration using Gaussian Quadrature 
[Ref. 6]. Thus, for example, the stiffness matrix is calculated by the 
triple summation: 


i i i 



where the weight coefficients and abscissa are given in Table 2. Note 
that the triple product matrix operation in Eq. (16) as well as the determinent 
of the Jacobian, |j|, must be evaluated at each integration point. This process 
could be very time consuming, and requires that efficient programming practices 
and mathematical techniques be used to minimize this time. In NASTRAN, for iso- 
tropic materials, the triple product in Eq. (21). is explicitly evaluated to 
avoid calculating zeros and thereby minimize the number of mathematical opera- 
tions performed. 


The number of integration points needed to evaluate the . stiffness , mass and 
load matrices depends on the element geometry, displacement function, and ma- 
terial property variations. Elements which are extremely distorted from a 
rectangular shape require more integration points. Best results, however, are 
obtained using rectangular elements as far as possible, and therefore extremely 
distorted elements should be avoided [Refs. 8 & 9]. It has been found that 
for most problems, satisfactory results may be obtained using a 2x2x2 
integration mesh for the linear element and a 3x3x3 integration mesh for 
the quadratic and cubic elements. These meshes are used as default values in 
NASTRAN. However, since good results have also been reported using smaller 
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■ t \/ j ■ . : 
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0.55555555555 

0.0 

0.88888888888 

±0.86113631159 

0.34785484514 

±0.33998104358 

0.65214515486 
















meshes than suggested above [Refs. 7 & 8], the user has the option of speci- 
fying the integration mesh size. 


,5.13.6; Transformation from Basic to Global Coordinates 

As previously stated, all computations for the isoparametric elements are 
carried out in the basic coordinate system. If the global coordinate system 
at any grid point is different from the basic system, NASTRAN transforms the 
final matrices and vectors into the global coordinate system using the appro- 
priate transformation matrix [T^] , corresponding to grid point i. This calcu- 
lation is identical to that performed for the other elements in NASTRAN. 

:5.13.7' Stress Recovery 

Element stresses may be obtained by combining Eqs. (4), (8) and (9) to yield: 

H = K] O] f^eH*e} [N] [T e } (22) 

where the matrices [C] and [N] are functions of the element coordinates £, 
r|, and £. in NASTRAN these stresses are calculated in the basic coordinate 
system and they are printed at the following locations, depending on element 
type: 

a. Linear Element - Eight corner points and at center of element. 

b. Quadratic and Cubic Element - Eight corner points, center of each 
edge, and at center of element. 

The principal stresses, principal angles, mean stress, and octohedral shear 
stress are also computed, and output at every point at which the basic stresses 
are computed. The mean stress, or hydrostatic pressure, is given by: 

o-n* -3 K + • (23) 
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The octahedral shear stress is given by: 



( 24 ) 


where S , S , and S are the three principal stresses, 
x y z 
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DIFFERENTIAL STIFFNESS AND BUCKLING 


7.4 ISOPARAMETRIC HEXAHEDRON THREE-DIMENSIONAL SOLID ELEMENTS 


The differential stiffness matrix for the isoparametric solid elements are. 
obtained by adding the energy due to initial stresses to the potential 
energy function. This additional energy is derived in Section 7.1, and is 
given by: 

'V = X [ [«s C«v + °>) +■ s (o * + "P + (a * + 

- A ^ ojy ~ ~ lUj } ] (LV 

where the rotations are given by the relations: 
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and may be expressed in terms of the grid point displacements by using Eq. (1) 
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Substituting Eq. (3) into (1) and adding this function to the Potential 
Energy expression yields the differential stiffness matrix: 


where 



f [ £ f i>ij 




'S 1 '- 



~ v i>t 


' x n 

° r *+ 


( 4 ) 


(5) 


As in the structural stiffness matrix, the evaluation of the integral in 
Eq. (4) is obtained by application of the Gaussian Quadrature Formula 
(Table 2, Section 5.13:5). 
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HEAT TRANSFER ANALYSIS 


8.2 VOLUME HEAT CONDUCTION ELEMENTS 

The volume hont conduction elements are the same as NASTRAN structural elements. The elements 
for which heat conduction is available are listed in the following table: 


Heat Conduction Elements 

Type 

Elements 

Linear 

BAR, R0D, C0NR0D, TUBE 

Planar 

TRMEM, TRIA1 , TRIA2, 
QDMEM, QUAD! , QUAD2 

Solid of Revolution 

TRIARG, TRAPRG 

Solid 

TETRA, WEDGE, HEXA1 , HEXA2 
IHEXI, IHEX2, IHEX3 


Scalar elements, single point constraints, and multipoint constraints are also available for heat 
transfer analysis. The same connection and property cards are used for heat transfer and struc- 
tural analysis. Linear elements have a constant cross-sectional area. For the planar elements, 
the heat conduction thickness is the membrane thickness. Elements with bending properties, such 
as BAR and TRIA1, have been included so that the user may use the same elements for the thermal 
and structural analyses of a given structure. The bending characteristics of the elements do not 
enter into heat conduction problems. The trapezoidal solid of revolution element, TRPRG, has been 
generalized to accept general quadrilateral rings (i.e., the top and bottom need not be perpendi- 
cular to the z-axis) for heat conduction only. 

The heat conduction elements are composed of constant gradient lines, triangles and tetra- 
hedra. The quadrilaterals are composed of overlapping triangles, and the wedges and hexahedra are 
formed from sub-tetrahedra in exactly the same way as for the structural case. The IHEXi 
elements are isoparametric hexahedron elements and are similar to the 
isoparametric solid elements described for structural analyses. 

Thermal conductivity and capacity are specified on MAT4 (isotropic) and MAT5 (anisotropic) 
bulk data cards. 

The heat conduction matrix for a volume heat conduction element may be derived from a thermal 
potential function in the same way that the stiffness matrix of a structural element Is derived 
from the strain energy function. The thermal potential function is 
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8.2.2 Three-Dimensional Isoparametric Solid Element Heat Transfer Matrices 


(add this new section) 

The heat transfer conduction matrix for the three-dimensional isoparametric 
elements [Refs. 1, 2, 3] are derived by using Equation (8). For these three 
elements, the temperature u at an exterior point is given by: 

u = [L ] {u } (23) 

e e 


where is a function of the curvilinear coordinates £, ri and £ shown in 
Figure 3, and are identical to the NL functions given in the structural 
analysis section. The derivatives of L g with respect to the basic Cartesian 
coordinates X, Y, Z are calculated in exactly the same manner as presented 
in the structural analysis section for these elements, i.e., 
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(24) 


where [J] is the Jacobian Matrix, and the derivatives of L_^ with respect to 

£, d and £ are listed in Tables 1 and 2 of the structural analysis section. 

Since, in general, the matrix [L . ] is a function of £, r| and £, the inte- 

® i 

gration of Equation (8) is carried out by Gaussian Quadrature numerical 
integration, and the heat conduction matrix is calculated from the equation: 


|j; 


(25) 
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where 



*L, 5L , 

ax ay az 

3L i * L i 

ax ay az 


( 26 ) 


' l n&p aL fj&p aL N&p 

ax ay az 


and NGP is the number of grid points in the element. NGP equals eight, 
twenty and thirty-two for the linear, quadratic and cubic isoparametric 
elements, respectively. The weight coefficients and abscissa are 
given in Table 1 in the structural analysis Section 5.13. 

Elements of the heat capacity matrix [B] are calculated by the coupled 
mass method, (see Section 5.13.4). The equation for this matrix is: 

[B] = / [L ] T [L] c dV (27) 

J V e e P 


where c^ is the heat capacity per unit volume. Examination of this equation 
shows that [B] is identical to the structural mass matrix with the single 
exception that the heat capacity c^ is used instead of the material density. 
The heat capacity matrix is also evaluated using the Gaussian Quadrature 
formula: 


[B] 


xxz 

*=• i'-' k--\ 
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(28) 


42 



REFERENCES 


1. 0. C. Zienkiewicz and Y. K. Cheung, The Finite Element Method in Struc- 
tura l and Continuum Mechanics , McGraw Hill Publishing Company, New 
York, 1967. 

2. 0. C. Zienkiewicz and C. J. Parekh, "Transient Field Problems: Two- 

Dimensional and Three-Dimensional Analysis by Isoparametric Finite 
Elements," International Journal for Numerical Methods in Engineering , 
January-March 1970. 

3. 0. C. Zienkiewicz and Y. K. Cheung, "Finite Elements in the Solution of 
Field Problems," The Engineer , pp. 507-10, September 24, 1965. 


43 






APPENDIX B 

USER'S MANUAL UPDATES 


The following sections contained in this appendix should be 
added to the NASTRAN User's Manual in thei r entirety: 

1.3.10 Isoparametric Solid Hexahedron Elements 
(including Figure 13) 

2.4.2 Bulk Data Card Descriptions: Add CIHEX1, 

CIHEX2 , CIHEX3, PIHEX, and PL0AD3 

The remaining sections are marked to indicate changes to the 
existing UM document. 
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The following new section and figure should be added to the User's Manual: 

I 

1.3.10 Isoparametric . Solid . Hexahedron Elements 

Three types of isoparametric solid hexahedron elements are provided for 
general solid structures. These elements (see Figure 13) are a linear, a 
quadratic, and a cubic isoparametric hexahedron. The theory is given in 
section 5.13 of the Theoretical Manual . These elements can be used with 
all other NASTRAN elements ,. except the axisymetric. elements . Connections 
are made only to the translational degrees of freedom at the grid points. 

The elements are defined by CIHEXI, CIHEX2, and CIHEX3 connection cards. 

All three of these cards reference the PIHEX property card. 

The isoparametric solid hexahedron elements are sophisticated elements. 
They allow the user to accurately define a structure with fewer elements 
and grid points than might otherwise be necessary with simple constant 
strain solid elements. The linear element generally gives best results 
for problems involving mostly shear deformations, and the higher order 
elements give good results for problems involving both shearing and bending 
deformations. Only a coupled mass matrix is generated to retain the inherent 
accuracy of the elements. Temperature, temperature dependent material proper- 
ties, displacements, and stresses may vary through the volume of the elements. 
The values at interior points of the element are interpolated using the 
isoparametric shape function. For best results, the applied grid point 
temperatures should not have more than a "gentle" quadratic variation in 
each of the three dimensions of the element. If the element has non-uniform 
applied temperatures, or if it is not a rectangular parallelopiped , 3 or 
more integration points should be specified on the PIHEX card. Severely 
distorted element shapes should be avoided. 

Stiffness, mass, differential stiffness, structural damping, conductance, 
and capacitance matrices may be generated with these elements. Piecewise 
linear analysis has not been implemented. 

The output stresses are given in the basic coordinate system. The 
stresses are assumed to vary through the element. Therefore, stresses are 
computed at the center and at each corner grid point of these elements. For 
the quadratic and cubic elements, they are also computed at the midpoint of 
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each edge of the element. In addition to the six normal and shear stresses, 
output. also includes . the principal stresses, the direction cosines of the 
principal planes, the mean stress . 


a 

n 


- t (a + a + a; ) 

3 x y z 


and the octahedral shear stress 
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where S 


x 


and are the three principal stresses. 


1.5.1 Static Loads 

The following pages should be modified as indicated. 


The minus sign is used to be consistent with the hydrostatic pressure 
computations for elements TETRA, WEDGE, HEXA1 and NEXA2. 
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STRUCTURAL MODELING 

Pressure loads on triangular and quadrilateral elements are defined with a PL0AD2 card. 

The positive direction of the loading is determined by the order of the grid points on the 
element connection card, using the right hand rule. The magnitude and direction of the load is 
automatically computed from the value of the pressure and the coordinates of the connected grid 
points. The load Is applied to the connected grid points. The PL0AD card Is used In a similar 
fashion to define the loading of any three or four grid points regardless of whether they are 
connected with two-dimensional elements. The PRESAX card is used to define a pressure loading on 
a conical shell element. 

Pressure loads on the isoparametric solid elements are defined with the 
PL0AD3 card. The pressure is defined positive outward from the element. 

The magnitude and direction of the equivalent grid point forces is automatic- 
' ally computed using the isoparametric . shape functions of the element to which 
the load has been applied. 

' ' . _ / ' i 

The GRAV card is used to specify a gravity load by providing the components of the gravity 
vector in any defined coordinate system. The gravity load is obtained from the gravity vector and 
the mass matrix assembled by the Structural Matrix Assembler (see Section 4.28 of the Programmer's 
Manual). The gravitational acceleration is not calculated at scalar points. The user is required 
to Introduce gravity loads at scalar points directly. 

The RF0RCE card Is used to define a static loading condition due to a centrifugal force 
field. A centrifugal force load Is specified by the designation of a grid point that lies on 
the'axls of rotation and by the components of rotational velocity in any defined coordinate 
system. In the calculation of the centrifugal force, the mass matrix is regarded as pertaining 
to a set of distinct rigid bodies connected to grid points. Deviations from this viewpoint, 
such as the use of scalar points or the use of mass coupling between grid points, can result In 
errors . 

Temperatures may be specified for selected elements. The temperatures for a R0D, BAR, C0NR0D 
or TUBE element are specified on the TEMPRB data card. This card specifies the average tempera- 
ture on both ends and. In the case of the BAR element. Is used to define temperature gradients 
over the cross section. Temperatures for two dimensional plate and membrane elements are specified 
on a TEMPP1 , TEMPP2, or TEMPP3 data card. The user defined average temperature over the volume Is 
used to produce in-plane loads and stresses. Thermal gradients over the depth of the bending 
elements, or the resulting moments, may be used to produce bending loads and stresses. 

( 

1.5-2 (7/1/70) Updated 
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STRUCTURAL MODELING 

element arc temperature-dependent by use of the MATTi card, they are always calculated from the 
"average" temperature of the element, i except for the isoparametric solid elements. 

The temperature is allowed to vary throughout the isoparametric solid elements. 

Hence, the material coefficients and the integration of the thermal loads require 

that the temperature at any point in the element be interpolated from the grid 

point temperatures specified. The mere presence of a thermal field does not imply the 
application of a thermal load. A thermal load will not be applied unless the user makes a specific 
request in the Case Control Deck. 

Enforced axial deformations can be applied to rod and bar elements. They are useful in the 
simulation of misfit and misalignment in engineering structures. As in the case of thermal expan- 
sion, the equivalent loads are calculated by separate .subroutines for each type of structural 

1 


\ 


1. 5-2a (7/1/70) Updated 
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2.4 Bulk Data Deck 

. N . / 

Section 2.4.2, Bulk Data Card Descriptions. Five new descriptions for the 
new bulk data cards should be inserted in this section. The new descriptions 
are contained in Figures 1 to 5. In addition, remark 4 for the TEMP, and 
TEMPD cards should be changed to read as follows: 

4. If the element material is temperature dependent, its 
properties are evaluated at the average temperature. 

The isoparametric solid elements are exceptions to this. 

Their properties are evaluated at the temperature com- 
puted by interpolating the grid point temperatures. 
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Input Data Card CIHEX1 - Linear Isoparametric Hexahedron Element Connection 

Description : Defines a linear isoparametric hexahedron element of the 

structural model. 

r 

Format and Example : 


1 2 3 4 5 6 7 8 9 10 


CIHEX1 

EID 

PID 

Gl 

G2 

G3 

G4 

G5 

G6 

abc 

CIHEX1 

137 

5 

3 

8 

5 

4 . 

9 

14 

ABC 


+bc 

G7 

G8 








+BC 

11 

10 









Field 


Contents 


EID Element identification number (Integer > 0) 

PID Identification number of a PIHEX property card (Integer > 0) 

G1,...,G8 Grid point identification numbers of connection points 



Remarks : 1. Element identification numbers must be unique with respect -to 

all other element identification numbers. 

2. Grid points Gl, G2, G3, G4 must be given in counter-clockwise 
order about one quadrilateral face when viewed from inside the 
element. G5, G6, G7, G8 are in order in the same direction 
around the opposite quadrilateral, with Gl and G5 along the 
same edge. 

3. There is no non-structural mass. 


4. The quadrilateral faces need not be planar. 

5. Stresses are given in the basic coordinate system. 

6. For structural problems, material must be defined by a MAT1 card. 

7. For heat transfer problems, material may be defined with either 
a MAT 4 or MAT 5 card. 


FIGURE 1 
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Input Data Card CIHEX2 - Quadratic Isoparametric Hexahedron Element Connection 


Description : Defines a quadratic isoparametric hexahedron element of the 

structural model. 

Format and Example: 


1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

CIHEX2 

EID 

PID 

G1 

G2~ 

G3 

G4 

G5 

G6 

abc 

CIHEX2 

110 

7 

3 

8 

12 

13 

14 

9 

ABC 


+bc 

G7 

G8 

G9 

G10 

Gil 

G12 

G13 

G14 

def 

+BC 

• 5 

4 

mm 

19 

20 

17 

23 

27 

DEF 


+ef 

G15 

G16 

G17 

G18 

G19 

G20 




+EF 

31 

32 

33 

28 

25 

24 





Field 


Contents 


EID Element identification number (Integer > 0) 

PID Identification number of a PIHEX property card (Integer > 0) 

G1,...,G20 Grid point identification numbers of connection points 
(Integer > 0, G1 £ G2 / ... ^ G20) 



FIGURE 2 
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Remarks 


: 1. Element identification numbers must be unique with respect to 

all other element identification numbers. 

2. Grid points G1,:..,G8 must be given in counter-clockwise order 
about one quadrilateral face when viewed from inside the element. 
G9,...,G12 and G13,...,G20 are in the same direction with Gl, G9 
and G13 along the same edge. 

3. There is no nonstructural mass. 

4. The quadrilateral faces need not be planar. 

5. Stresses are given in the basic coordinate system. 

6. For structural problems, material must be defined by 
a MAT1 card. 

7. For heat transfer problems, material may be defined 
with either a MAT 4 or MAT5 card. 


FIGURE 2 (cont'd) 



Input Data Card CIHEX3 - Cubic Isoparametric Hexahedron Element Connection 

Description : Defines a cubic isoparametric hexahedron element of the 

structural model. 

Format and Example : 


1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

CIHEX3 

EID 

PID 

Gl 

G2 

G3 

G4 

G5 

G6 

abc 

CIHEX3 

15 

3 

4 

9 

12 

17 

18 

19 

ABC 


+bc 

G7 

G8 

G9 

G10 

Gil 

G12 

G13 

G14 

def 

+BC 

20 

13 

10 

7 

6 

5 

22 

25 

DEF 


4-ef 

G15 

G16 

G17 

G18 

G19 

G20 

G21 

G22 


+EF 

26 

23 

28 

31 

32 

29 

36 

41 

GHI 


+hi 

G23 

G24 

G25 

G26 

G27 

G28 

G29 

G30 

jkl 

+HI 

44 

49 

50 

51 

52 

45 

42 

39 

JKL 


+kl 

G31 

G32 








+KL 

38 

37 









Field 


Contents 


EID 

PID 

Gl, . . . ,G32 


Element identification number (Integer > 0) 

Identification number of a PIHEX property card (Integer > 0) 
Grid point identification number of connection points 
(Integer >0', Gl ^ G2 ^ ... / G32) 


FIGURE 3 
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Remarks : 1. Element identification numbers must be unique with respect to 

all other element identification numbers.' 

2. Grid points G1,...,G12 must be given in counter-clockwise order 
about one quadrilateral face when viewed from inside the element. 
G13,...,G16; G17,...,G20; and G21,...,G32 are in the same direc- 
tion with Gl, G13, G17, G21 along the same edge. 

3. There is no nonstructural mass. 

4. The quadrilateral face need not be planar. 

5. Stresses are given in the basic coordinate system. 

6. For structural problems, material must be defined by 
a MAT1 card. 

7. For heat transfer problems, material may be defined 
with either a MAT4 or MAT5 card. 


FIGURE 3 (cont'd) 
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Input Data Card PIHEX - Isoparametric Hexahedron Properties 


Description : Defines the properties -of an isoparametric solid element. of the 

structural mode, including. a material reference and the number 
of integration points. Referenced by the CIHEX1, CIHEX2, and 
CIHEX3 cards. 

Format and Example: . 


1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

PIHEX 

PID 

MID 

CID 

NIP 

AR 

ALFA 

BETA 



PIHEX 

15 

3 


3 



5.0 




Field 


Contents 


PID 

MID 

CID 

NIP 

AR 

ALFA 

BETA 


Property identification number (Integer > 0) 

Material identification number (Integer > 0) 

Identification number of the coordinate system in which the 
material referenced by MID is defined (Integer > 0 or blank) 

Number of integration points along each edge of the element 
(Integer = 2, 3, 4 or blank) 

Maximum aspect ratio, (ratio of longest to shortest edge) of the 
element (Real > 1.0 or blank) 

Maximum angle between the normals of two sub-triangles comprising 
a quadrilateral face (Real, 0.0 < ALFA < 180.0, or blank) 

Maximum angle between the vector connecting a corner point to an 
adjacent midside point and the vector connecting that midside 
point and the other midside or corner point (Real, 0.0 < BETA^_ 180.0, 
or blank) 


Examples of Field Definitions : 



FIGURE 4 
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Remarks : 


. All PIHEX cards must have unique identification numbers. 

. CID is not used for isotropic materials. 

. The default for CID is the basic coordinate system. 

. The default for NIP is 2 for IHEXI and 3. for IHEX2 and IHEX3. 

. AR, ALFA, and BETA are used for checking the geometry of the 
element. The defaults are: 



AR 

ALFA 

(degrees) 

BETA 

(degrees) 

C IHEXI 

5.0 

45.0 

- 

CIHEX2 

10.0 

45.0 

45.0 

CIHEX3 

15.0 

45.0 

45.0 


FIGURE 4 (cont'd) 
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Input Data Card PL0AD3 - Pressure Load on a Face of an Isoparametric Element 

Description : Defines a uniform static pressure load applied to a surface of 

an isoparametric element only. 

Format and Example : 


1 

2 

3 

r T~ 

5 

6 7 

8 

A 

9 

10 

PL0AD3 

SID 

P 

EIDl 

Gil 

G12 

EID2 

G21 

G22 


PL0AD3 

3 

-15.1 

15 

7 

25 

16 

117 

135 



Field 


Contents 


SID 

P 

EID1 V 
EID2 / 

Gil, G12 ^ 
G21, G22 / 


Load set identification number (Integer > 0) 
Pressure value (Real, force per unit area) 

Element identification number (Integer > 0) 


Grid point identification number of two grid points at diagonally 
opposite corners of the face on which the pressure acts 
(Integers > 0) 


Remarks : 


1. Load sets must be selected in the Case Control Deck (LOAD = SID) 
to be used by NASTRAN. 

2. At least one EID must be present on each PL0AD3 card. 

3. All elements referenced must exist. 

4. Computations consider the pressure to act positive outward on 
specified face of element. 


FIGURE 5 




4.2 Structure Plotting 


Section 4. 2. 2.1 SET Definition Cards, the following element types are also 
now permissible: 

IHEX1, IHEX2, IHEX3, TETRA, HEXAl, HEXA2 , WEDGE 
Section 4. 2. 2. 3 PL0T Execution Card, the plot labels for the solid elements are 


Element Type 

Label 

IHEX1 

XL 

IHEX2 

XQ 

IHEX3 

XC 

TETRA 

TT 

HEXAl 

HI 

HEXA2 

H2 

WEDGE 

WG 
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6.2 NASTRAN System. and User Messages 

Section 6.2.3, Functional Module Messages, the following new messages should 
be added to this section: 

2141 *** USER FATAL MESSAGE 2141, IHEX* ELEMENT . NUMBER. *** INSUFFICIENT 

C0RE T0 C0MPUTE ELEMENT MATRIX 

2142 *** USER FATAL MESSAGE 2142, IHEX* ELEMENT NUMBER.*** ILLEGAL 

GE0METRY , text 

The type of geometry error is identified in "text". The possibilities 
are : 


AR: EXCEEDED \ 
ALFA EXCEEDED > 
BETA EXCEEDED ) 
REVERSED NUMBERING 


Either correct the. element or increase the 
allowable value on the PIHEX card for this 
element . 


The element was numbered in a clockwise fashion 
rather than counter-clockwise as required. This 
would result in a left-handed element coordinate 
system. Correct the numbering sequence on the 
CIHEXi card for this element. 


C00RDINATES 0F TW0 The coordinates of all connections of the element 
P0INTS ARE THE SAME must be different. 


2143 *** USER FATAL MESSAGE 2143, SINGULAR JAC0BIAN MATRIX F0R IS0PARAMETRIC 
ELEMENT NUMBER *** 


The element is severely warped. or the outer surface: of the element is 
folded through itself. Check the connection card for this element and 
the coordinates of the points it connects. 

2144 *** USER FATAL MESSAGE . 2144 , PL0AD3 CARD FR0M L0AD SET *** REFERENCES 
MISSING 0R N0N-IS0PARAMETRIC ELEMENT *** 

2145 *** USER FATAL MESSAGE .2145, INVALID GRID P0INTS 0N PL0AD3 CARD F0R 
ELEMENT *** 

Either the element does not connect the specified grid points, or the grid 
points do not identify the diagonal of a face of the element. 

4024 *** USER FATAL MESSAGE 4024, STRESSES REQUESTED F0R SET *** WHICH 
C0NTAINS N0 VALID ELEMENT ID-S 

The set of elements for which stresses were requested in this subcase 
contains only ID's for nonexistent elements. 
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2. DATA BLOCK AND TABLE DESCRIPTIONS 


2.3.2 Data Blocks Output from Module IFP 


Section 2. 3. 2. 2 GE0M2 (TABLE), the following data should be added to this 
section: 

I Card Types and Header Information : 


Header Word 1 Header Word 2 Header Word 3 

Card Type Trailer Bit Position Internal Card Number 


CIHEX1 7108 


71 

251 

CIHEX2 7208 


72 

252 

CIHEX3 7308 


73 

253 

Card Type Formats : 

CIHEX1 (10 words) 

EID 

PID 

G1 


G2 

. . . 

G8 

CIHEX2 (22 words) 

EID 

PID 

G1 


G2 

• • • 

G20 

CIHEX3 (34 words) 

EID 

PID 

G1 


G2 

• • • 

G32 

Section 2. 3. 2. 3 GE0M3 (TABLE) 

the 

following data 

should be added to this 

section: 

| Card Types and Header Information: 


Header Word 

1 

Header Word 2 

Header Word 3 


Card Type Card Type Trailer Bit Position Internal Card Number 


PL0AD3 7109 

Card Type Formats : 

PL0AD3 (5 words) SID 

G1 


71 


P EID 

G2 


255 


Section 2. 3. 2. 5 EPT (TABLE), the following data should be added to this 
section: 


Card Types and Header Information : 

Header Word 1 Header Word 2 Header Word 3 

Card Type Card Type Trailer Bit Position Internal Card Number 


PIHEX 7002 
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Card Type Formats: 




PIHEX (7 words) 

PID 

MID 

CID 


NIP 

AR 

ALFA 


BETA 

' 65 ' 


2.3.7 Data Blocks Output from Module GP3 


Section 2. 3. 7.1 SLT (TABLE), the following data should be added to this 
section under "Notes", item 2: 


13 = PL0AD3 



Word 

T yp e 

Item 

1-6 

R 

Pressures 

7-38 

I 

Internal grid numbers of which the 
last 12 or the last 24 may be zero 

Under "Notes", item 

3 should begin: 


3. With the exception of GRAV, PL0AD, and PL0AD3 card types, 


2.3.8 Data Blocks Output from Module TA1 
The EST table should be amended as follows: 
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2.3.8. 1 EST (TABLE) . 

Description 

-Element Summary Table. 

The EST is a collection of data for all elements of the structural model. It contains one 
logical record for each element type. For each element: connection data, properties data, 


basic grid point data and the element temperature-sare grouped, 
that belong to super elements are not Included in /the EST. 

T . . c , cLa£a^ 

Table Format 


General elements and elements 


Record 

0 

1 


Word 


Type 


I 


1 

2 - 1+1 
i+2-i+.i+l 
i+j+2-i+j+k+l . 

i+j+k+2-JL Je + rn + 1 


n+1 


Header record 

Element type 
ECT section 
EPT section 
BGPDT section 
Elcm e nt - tompoP d t u ne 
ETT 

End-of-file 


Item 


repeated 

for 

each 

element 


repeated 
for 
each 
el ement 
type 


Notes 


1. i = number of words in ECT section, 

j = number of words in EPT section. 

It = number of words in BGPDT section. 

'A - -rvumJbtA. wVLit iAt- ETTaajCJ>>k.. 

2. The number of records in the EST corresponds to the number of separate element types In 
the model . 


Summary of EST Formats 


3. The EST is generated in subroutine TA1A. 

ETT ctum. 

ECT Section EPT Section BGPDT Section TUo'ulA, Total 

ErieweM' Words 

number Number 


Element 


numuer 


Per 


Jyp.e 

Mnemonic 

of Words 

of Words 

of Words 

-atwe 

Element 

1 

RpD 

3 

5 

8 

Yos 

17 

2 

BEAM 

19 

19 

8 

Yis 

47 

3 

TUBE 

3 

4 

8 

Y|s 

16 

4 

SHEAR 

5 

3 

16 

Yes 

25 

5 

TWIST 

5 

3 

16 

Yas 

25 

6 

TRIA1 

5 

9 

12 


27 

7 

TRBSC 

5 

7 

12 

Yfis 

25 

8 

TRPLT ■ 

5 

7 

12 

Yis 

25 

9 

TRMEM 

5 

3 

12 

Yfo 

21 

10 

C(?NR(?D 

8 

0 

8 

Ybs 

17 

11 

ELAS1 . 

5 

3 

0 

t;d 

8 

12 

ELAS2 

8 

0 

0 

N$ 

8 

13 

ELAS3 

3 

3 

0 

No 

6 

14 

EIAS4 

4 

0 

0 

fja 

4 

15 

flCPLT 

6 

7 

16 

yes 

30 

16 

QDMEM 

6 

3 

16 

Y«S 

26 

17 

TRIA2 

5 

3 

12 

Yis 

21 


C 




Yes t, 

Wo "to o 
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ECT 

section 

for element 

type = 62: 





Word 

T yP e 

Item 





1 

I 

Element ID 





2-9 

I 

SIL numbers for 

grid 

points 

1-8 

ECT 

section 

for element 

type = 63: 





Word 

Type 

Item 





1 

I 

Element ID 





2-21 

I 

SIL numbers for 

grid 

points 

1-20 

ECT 

section 

for element 

type = 64: 





Word 

Typ_e 

Item 





1 

I 

Element ID 





2-33 

I 

SIL numbers for 

grid 

points 

1-32 
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2.3.8 Data Blocks Output From Module TA1 (cont'd) 


EPT 

section 

for element 

type = 62, 63, 64: 




Word 


Item 




1 

I 

Material ID 




2 

I 

CID 




3 

I 

NIP 




4 

R 

AR 




5 

R 

ALFA 




6 

R 

BETA 



ETT 

section 

for element 

type = 1-61: 




Word 

Type 

Item 




1 

R 

Average element temperature 

ETT 

section 

for element 

type = 62: 




Word 

T yP e 

Item 




1-8 

R 

Temperature at grid 

points 

1-8 

ETT 

section 

for element 

type = 63: 




Word 

Type 

Item 




1-20 

R 

Temperature at grid 

points 

1-20 

ETT 

section 

for element 

type = 64: 




Word 

T yP e 

Item 




1-32 

R 

Temperature at grid 

points 

1-32 
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2.3.51 Element Stress Output Data Description 
The following data should be added to this section: 


Element 
Type Name 
61 CIHEX1 + 


62 CIHEX2 + 

63 CIHEX3 t 


Real Element Stresses Complex Element Stresses 

Word. or Word or Real 

Component Item Component Item Imag. 


2 

External grid point ID 

2 

External grid 





point ID 


3 

Normal - x 

3 

Normal - x 

R 

4 

Shear - xy 

4 

Normal - y 

R 

5 

First principal 

5 

Normal - z 

R 

6 

First principal x cosine 

6 

Shear - xy 

R 

7 

Second principal x cosine 

7 

Shear - yz 

R 

8 

Third principal x cosine 

8 

Shear - zx 

R 

9 

Mean stress 

9 

Normal - x 

I 

10 

Octahedral shear stress 

10 

Normal - y 

I 

11 

Normal - y 

11 

Normal - z 

I 

12 

Shear - yz 

12 

Shear - xy 

I 

13 

Second principal 

13 

Shear - yz 

I 

14 

First principal y cosine 

14 

Shear - zx 

I 

15 

Second principal y cosine 




16 

Third principal y cosine 




17 

Normal - z 




18 

Shear - zx 




19 

Third principal 




20 

First principal z cosine 




21 

Second principal z cosine 




22 

Third principal z cosine 





Note CIHEX1 


Note CIHEX1 


2 

First external grid 

2 

First external 



point ID 


grid point ID 


3 

Normal - x 

3 

Normal - x 

R 

4 

Shear - xy 

4 

Normal - y 

R 

5 

First principal 

5 

Normal - z 

R 

6 

First principal x cosine 

6 

Shear - xy 

R 

7 

Second principal x cosine 

7 

Shear - yz 

R 

8 

Third principal x cosine 

8 

Shear - zx 

R 

9 

Mean stress 

9 

Second external 





grid point ID 


10 

Octahedral shear stress 

10 

Normal - x 

I 

11 

Second external grid 

11 

Normal - y 

I 


point ID 




12 

Normal - y 

12 

Normal - z 

I 

13 

Shear - yz 

13 

Shear - xy 

I 

14 

Second principal 

14 

Shear - yz 

I 

15 

First principal y cosine 

15 

Shear - zx 

I 

16 

Second principal y cosine 




17 

Third principal y cosine 
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2.3.51 Element Stress Output Data Description (cont'd) 


Element Real Element Stresses 

Word or 

Type Name Component Item 


Complex Element Stresses 
Word or Real 

Component Item Imag. 


18 Normal - z 

19 Shear - zx 

20 Third principal 

21 First principal z cosine 

22 Second principal z cosine 

23 Third principal z cosine 


The stresses are repeated for each stress point within each element. 
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Element 

Force Output 

Data Descripti 

on 



following data should be 

added to this 

section: 



Element 

Real Element Stresses 

Complex 

Element Stresses 



Word or 


Word or 


Real 

Type 

Name 

Component 

Item 

Component 

Item 

Imag 

62 

CIHEX1 

- 

Undefined 

- 

Undefined 

- 

63 

CIHEX2 

- 

Undefined 

- 

Undefined 

- 

64 

CIHEX3 

— 

Undefined 

— 

Undefined 

_ 


2.4.2 Executive Tables Not Permanently Core Resident 

Section 2. 4. 2. 8 IFPX1 (Master Card Name Table), the following data should 
be added to this section: 


Word No. 
In IFPX1 

Bit No. 
In IFPX0 

Contents 

Output 
Only (PMR) 

Supported 

502 

251 

CIHEX1 

- 

- 

504 

252 

CIHEX2 

- 

- 

506 

253 

CIHEX3 

- 

- 

508 

254 

PIHEX 

- 

- 

510 

255 

PL0AD3 

- 

- 
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3. SUBROUTINE DESCRIPTIONS 


3.4.36 PREMAT (Material Property Utility) 

Section 3.4.36.3 Calling Sequence, common block /MATOUT/ has been lengthened 
to 21 words. Its description should be altered as follows: 

C0MM0N/MAT0UT/ - (Output Common Block). Length 21 words. Depending 

upon the values of INFLAG, the output common block is 
defined variously as follows: 


1. MAT! Format (INFLAG = 1) 


Word 


Symbol 



8 

9 



10 o $ 

11-20 


21 TDEP 


Definition . , 

Young's modulus (modulus of elasticity) ^ 

f i 

Shear Modulus jj 

Poisson's ratio ! J 

Density 

Thermal expansion coefficient 
Thermal expansion reference temperature 

Structural element damping coefficient 

Stress limit for tension 

Stress limit for compression 

Stress limit for shear . j 

Undefined 1 

Temperature dependence flag (logical) 
T - material referenced was 
temperature dependent 
F - material referenced was not 
temperature dependent 


For the tables in items 2 through 8, corresponding to values of INFLAG from 
2 to 8, the entries like 

| 16-20 - Undefined 

should be changed to 

I 16-21 - Undefined 
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4. MODULE FUNCTIONAL DESCRIPTIONS 


4 . 5 Executive Preface Module IFP (Input File Processor) 


The following data should be added to Table 1 and Table 2: 


A 

B 

C 

D 

E 

F 

G 

H 

I 

J 

K 

L 

M 

N 

0 

IJHK 

251 

CIHEX1 

8 

GE0M2 

-2 

12 

16 

951 

1 

7108 

71 

S5 

5000 

-1 

33B5 

252 

CIHE-X2 

8 

GE0M2 

-2 

24 

28 

951 

1 

7208 

72 

S5 

5100 

-1 

33B6 

253 

CIHEX3 

8 

GE0M2 

-2 

36 

4 

951 

1 

7308 

73 

S5 

5200 

-1 

33C1 

254 

PIHEX 

2 

EPT 

-2 

4 

12 

981 

1 

7002 

70 

S5 

5300 

-1 

33C2 

255 

PL0AD3 

9 

GE0M3 

-2 

8 

12 

949 

0 

7109 

71 

S5 

5400 

-1 

33C3 


4.24 Functional Module PL0T (Structural Plotter) 

Section 4.24.8.8 Subroutine name: SHAPE, the last sentence of item 4 

should read: 

As each element is read, it is drawn, taking into account whether 
the element is one-, two-, or three-dimensional. 
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4.25 Functional Module GP3 (Geometry Processor-Phase 3) 

Section 4.25.8.2 Subroutine: GP3C, this entire section should be replaced 

with the following: 

1. Entry Point: GP3C 

2. Purpose: 

To convert PL0AD2 data (if present) to PLOAD format, merge PL0AD2 
data with PLOAD data (if present) and write the resulting data on 
SCR2, a scratch file. Then, to convert PL0AD3 data (if present) 
to an expanded format and write the resulting data on SCR2. 

3. Calling Sequence: CALL GP3C 

4 . Method : 

A. The trailer bits for PL0AD2 and PL0AD3 cards are checked. If 
no PL0AD2 cards exist, Step B is executed. Otherwise, PL0AD2 
cards are read into core from GE0M3. Six words are used for 
each entry. The first word (set identification) is set nega- 
tive and the sixth word of each entry is set to zero. Step C 
is executed. 

B. If no PL0AD3 cards exist* a return to GP3 is issued. Otherwise, 
PL0AD3 cards are read into core from GE0M3. 39 words are used 
for each entry, which contains all pressures on one element. 

As each PL0AD3 card is read, the element ID is checked against 
those of PL0AD3 cards previously read. If a match is found, the 
data is added to that entry (pressures on the same face are 
summed). If no match is found, the data is stored in the next 
unused entry. The first word (set identification) is set nega- 
tive and unused words are set to zero. 

C. GE0M2 is opened and the header record is skipped. The following 
steps occur for each record on GE0M2. 

Cl. The three-word header is read. /GPTA1/ (see Section 2.5, P.M.) 
is searched for a match. If no match is found, the record is 
skipped and the process is repeated. If an end-of-file is 
encountered, Step 4 is executed. If a match is found, at 
test on element type is made. If not a two-dimensional 
element and PL0AD2 data is being processed or if not an 
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isoparametric solid element and PL0AD3 data is being 
processed, the record is skipped and the process is 
repeated. Otherwise, step C2 is executed. 

C2. An entry of the current element type is read. A linear 
search through the PL0AD2 or PL0AD3 data is made to find 
a match on element identification. If no match is found, 
the next entry is read. For each match which is found 
for PL0AD2 data, the grid identification numbers which 
connect the element are stored in the corresponding 
PL0AD2 entry and the first word of that entry is set 
positive. For each match which is found for PL0AD3 data, 
the pressure factors are rearranged according to element 
face, the grid identification numbers which connect the 
element are stored in the corresponding PL0AD3 entry, and 
the first word of that entry is set positive. 

D. A pass through each entry in the PL0AD2 or PL0AD3 data is made. 

For each entry for which the first word is negative, an error 
message is queued or written and the NOGO flag is turned on. If 
upon completion of the pass, the NOGO flag is on and all PL0AD2 
and PL0AD3 data has been processed, PEXIT is called. If the NOGO 
flag is on and all data has not been processed, control is passed 
to Step B. Otherwise Step E is executed. 

E. If PL0AD3 data is being processed, Step F is executed. Otherwise 
LOCATE is called to position GE0M2 to PLOAD data. If none exists. 
Step F is executed. Otherwise, the PLOAD data is read into core 
following the PL0AD2 data. The combined list is sorted by SORT 
on set identification number. 

F. The data in core is written on SCR2 (All PL0AD2/PL0AD data is 
written as one logical record and all PL0AD3 is written as the 
next logical record). If all PL0AD3 data has been processed, a 
return to GP3 is given. 
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Allocation of open core storage in GP3C is as follows: 


1. For processing PL0AD2 data: 


NPLD2 


BUF1 


Before conversion 


PL0AD2 set ID 


Scale factor 


Element ID 


Same format as 1st 
entry 


GIN0 buffer 


C0MM0N/GP3C0R/Z(lj 


| Six words 
per 

PL0AD2 
I card 


Last 
1 PL0AD2 
entry 


NPLD2 


BUF1 


After conversion 


+ PL0AD2 or PL0AD set ID 


Scale factor 


Grid point ID-j 


Grid point ID 3 


Grid point ID 3 


Grid point ID^ or 0 


Same format as 1st entry 


GIN0 buffer 


Six words 
per entry 


\ 


, \ 


, / 


- > 


2. For processing PLOAD3 data: 
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C0MM0N/GP3C0R/Z(1) 


NPLD2 

BUFl 


Before Conversion 


After Conversion 


- PL0AD3 set ID 

Element ID 

Pressure 1 


Pressure 12 

Grid point ID^ 
for pressure 1 

Grid point ID^ 
for pressure 1 

• 

• 

• 

• 

• 

• 

Grid point ID^ 
for pressure 12 

Grid point ID£ 
for pressure 12 

0 

• 

Same format as 
first entry 

GIN 0 buffer 


39 words 
per entry 



39 words 
per entry 
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Section 4.25.8.3 Subroutine: GP3A, the parenthetical expression in item 4 

in the second sentence should read: 

| (or SCR2 for PL0AD and PL0AD3 data) 

Also for item 4, the parenthetical expression in the next-to-last paragraph 
should read : . 

I (except GRAV, PL0AD, and PL0AD3 data) 


Section 4.25.9.1 Allocation of Core Storage, the entry for GP3C should be 
changed to : 

GP3C: Maximum requirement = MAX [6* (number of PL0AD2 + number of 

PL0AD cards), 39* (number of elements with PL0AD3 pressures)] 

+ one GIN0 buffer. 


Section 4.25.10 Diagnostic Messages, the following messages may additionally 
be issued by GP3: 

1 2144, 2145. 


Figure 1 - Flowchart for module GP3, this figure should be modified as shown: 




















4 . 26 Functional Module TA1 (Table Assembler) 


Section 4. 26. 7.2 TA1A, this section should be replaced by: 

Assembly of the Element Summary Table is performed in two steps. 

For the first step, the EPT is read into core one property at a time. 

The ECT is read one element at a time. For each element the referenced 
property data are found by performing a binary search in the EPT in core. 
The ECT and EPT, data are written on SCR1, a scratch file, one element at 
a time, one logical record per element type. 

To initiate the second step, the BGPDT and SIL data blocks are read 
iin;to core. Data from SCR1 are read one element at a time. Internal 
indices for the grid points are used as pointers into the BGPDT and SIL 
tables. The temperature data for the element is extracted from the GPTT 
data block using subroutine TA1ETD. The internal indices are now replaced 
with corresponding scalar index values. A line comprising ECT, EPT, 

BGPDT and GPTT data for the element is written on the EST. Each logical 
record of the EST comprises all data of one element type. 
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Allocation of core storage during the second step is as follows: 

C0MM0N/TAA1/Z(1 ) 

. L0CBGP=1 


NBGP 

N8GP+1 


+WP- 


6UF3 

BUF2 


BUF1 


Coordinate system ID 

X 

y 

z 

• 

SILj 

SIL 2 


JLkmont jfcypo-.aad-N ■ 

vfi ym&nMS-- 


• 

• 

• 

GIN0 buffer C J- 

GIN0 buffer 

GIN0 buffer 


4 words/entry' 


} 1 word/entry 


-& H w or ds/eatry] 


BGPDT table 

one entry per point 

in problem 


SIL table 

one entry per point 
in problem 

-eae-e&t<yt-p 
defe i nad e4 ement 


Section 4.26.7.3 TA1B, the fifth sentence in the next to last paragraph 
should read: 

If the BGPDT and SIL can be held in core 
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The following section should be added: 

4.26.8.8 Subroutine Name: TA1ETD 

1. Entry Point: TA1ETD 

2. Purpose: To extract the temperature data for an element from 

the GPTT data block. 

3. Calling Sequence: CALL TA1ETD (ELID, TI, GRIDS) 

ELID - ID of the element for which temperature data is 
desired - integer-input 

TI - Array of temperature data - real-output 
.GRIDS - Number of grid points of the element. If GRIDS=0, 
.'average element temperature is returned in TI(1), (j 
integer-input ; 

Section 4.26.9.1 Allocation of Core Storage, step (2) for TA1A should read: 

Step (2) : Maximum core storage equals 5* (number of grid and scalor. 

points in model) plus three GIN0 buffers. 

Section 4.26.10 Diagnostic Messages, the following message may additionally 
be issued by TA1: 


4016 



4.27 Functional Module SMA1 (Structural Matrix Assembly - Phase 1) 

The following two sections should be; added: 

4.27.8.42 Subroutine Name: KIHEX 

1. Entry Point: KIHEX 

2. Purpose: To generate the element stiffness matrix for an 

IHEX1 , IHEX2, or IHEX3 element. 

3. Calling Sequence: CALL KIHEX (ITYPE) 

( 1 - IHEX1 j 

ITYPE = < 2 - IHEX2 > - integer-input 

I 3 - IHEX3 ) 


4.27.8.43 Subroutine Name: IHEXSD 

1. Entry Point: IHEXSD 


2. Purpose: To generate isoparametric shape function data. 


3. 


Calling Sequence: CALL IHEXSD (ITYPE, SHP, DSHP, JINV, 

DETJ, EID, XI, ETA, ZETA, C0RD>. 

/ / 


ITYPE 

SHP 

DSHP 


JINV 

DETJ 

EID 

X! | 
ETA > 
ZETA j 

C0RD 


- same as in KIHEX calling sequence - ,'integer-input; ) 

- array of values of the shape functions - double 
precision-output . 

- array of values of the derivatives with respect to 
element coordinates of the shape functions - double 
precision-output 

- array of the inverse of the Jacobian matrix - 
double precision-output 

- determinant of the Jacobian matrix - double 
precision-output 

- Element ID - integer-input 

Element coordinates at which shape functions are 
evaluated - double precision-input 

- array of basic coordinates of the grid points of 
the element - real-input 
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Section 4.27.10 Diagnostic Messages, the following additional user fatal 
error messages may be issued by the isoparametric solid element routines: 

I 2141, 2142, 2143. 


4 . 28 Functional Module SMA2 (Structural Matrix Assembler - Phase 2) 
The following two sections should be added: 

4.28.8.27 Subroutine Name: MIHEX 

1. Entry Point: MIHEX 

2. Purpose: To generate the coupled mass matrix 

for an IHEX1, IHEX2 , or IHEX3 element. 

3. Calling Sequence: CALL MIHEX (ITYPE) 

(l - IHEX1 ] 

ITYPE = < 2 - IHEX2 > - integer-input 
( 3 - IHEX3 ) 

4.28.28 Subroutine Name: IHEXSD 

See Section 4.27.8.43 
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4 . 41 Functional Module SSG1 (Static Solution Generator - Phase 1) 


Section 4.41.7 Overview of the Method Used in SSG1, the second sentence in 
item 1 should read: 

Pressure loads may be applied to an area defined by three or four 
grid points or to a face of an isoparametric solid element. 


Section 4.41.8 Direct Applied Loads, the first sentence of this section 
should read : 

Direct loads are applied to the structural model by means of 
F0RCE , F0RCE1, F0RCE2 , GRAV, M0MENT , M0MENT1, PL0AD , PL0AD2 , PL0AD3 , 
RF0RCE, and SL0AD Bulk Data Cards and the PRESAX card which is used 
for the axisymmetric conical shell problem only. 


The following four sections should be added: 

4.41.8.9 PL0AD3 Card .Processing . \ 

The data contents from a PL0AD3 card, as converted in GP3, are: 


P- ,P_ , . . . ,P, = Normal pressure value on faces 1 to 6 

go g 

i’ N 2’ ' ' * ,N 32 = I nterna l grid point numbers defining element on 

which pressure acts (N® to ma Y be zero depend- 
ing on element type) . 

The pres sure- To ad T vectbr T -on"the^i— element grid point is: 

6 2 2 

{f.} = 22 p k Z2 12 n® {b > 

k=l m=l n=l • 


til S 

Where P^ is the normal pressure on the k face, is the isoparametric shape 

function, and {B } is an area factor times a vector of direction cosines 
^ th 

defining the outward normal to the k face at each integration point. The 
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summation on m and n are for integration over the area of the k face using 

the method of Gaussian Quadrature. The load vector transformed to global 

th 

coordinates for the i point is: 


• {P } = [T/'] T { F . } 

8 ± i i 

where [T^] is the basic to global transformation matrix. 

The computational method used in computing the load vector is as follows: 

Each element of a vector { F } of length three times the number of grid points 
is set to zero. The following steps are repeated at each integration point 

for each face which has a non-zero applied pressure: 

/ 

s —I 

1. Isoparametric utility routines are called to compute N^'s and [J] 

(See Section 4.87.16.2). 

2. For each grid point i, and element face k, { B } is multiplied by N S P 

K. 1 rC 

and added to {F.}. (B } is simply +1 times det [J] times a column 
-1 lk 

of [J] . The column and sign are dependent on the face numbers as 

follows : 



to the global static load vector {P }. 
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4.41.11.38 Subroutine Name: PL0AD3 

1. Entry Point: PL0AD3 

2. Purpose: To apply loads due to PL0AD3 cards 

3. Calling Sequence: CALL PL0AD3 

C0MM0N/L0ADX/ 

L0ADX - See description of /L0ADX/ above (section 4.41.11.8) 

4.41.11.39 Subroutine Name: IHEX 


1. Entry Point: IHEX 

2. Purpose: To calculate an element thermal load vector for the 

IHEX1, IHEX2, and IHEX3 elements in the SSG1 module. 

3. Calling Sequence: CALL IHEX (TEMPS, PG, TYPE) 


C0MM0N/TRIMEX/ 


TEMPS - Array of grid point temperatures - real-input 
- Load vector array - real-input 


PG 

TYPE 

TRIMEX 


( 1 - IHEXIA . , 

<2 - IHEX2 > - integer-input 
(.3 - IHEX3 ) 

EST entry for the element - mixed-input 


4.41.11.40 Subroutine Name: QIHEX 


1. Entry Point: QIHEX 

2. Purpose: To calculate an element heat generation vector for 

the isoparametric solid elements for a unit external 
heat per unit volume. 

3. Calling Sequence: CALL QIHEX (TYPE,P) 

( 1 - IHEX1 ) 

TYPE = < 2 - IHEX2 V - integer-input 
( 3 - IHEX3 ) 

P - Array containing force vector - real-output 

Section 4.41.13 Diagnostic Messages, the following message may additionally 
be issued by SSG1: 

2143 
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4 . 46 Functional Module SDR2 (Stress Data Recovery - Phase 2) 


Section 4.46.8.21 Subroutine Name: SDR2D, paragraph 3 of item 4 should 

read: 

Core and GIN0 buffers are then allocated as required for a) the Case 
Control data block, b) the Element Deformation Table, c) the Grid 
Point Temperature Table, d) the element stress matrices, and e) the 
EQEXIN Data Block. If there is insufficient space in core for the 
element stress matrices, they are maintained on the scratch data 
block generated in stage III. 


The following three sections should be added. 

4.46.8.37 Subroutine Name: SIHEX1 

1. Entry Point: SIHEX1 

2. Purpose: To generate element stress matrices for the IHEX1, 

IHEX2, and IHEX3 elements, SIHEX1 must be called 
once for each stress point in an element. 

' 3. Calling Sequence: CALL SIHEX1 (TYPE, STRSPT) 

( 1 - IHEX1 } 

TYPE = < 2 - IHEX2 - integer-input 
( 3 - IHEX3 ) 

STRSPT - stress point number - integer-input 

4.46.8.38 Subroutine Name: IHEXSS 

1. Entry Point: IHEXSS 

2. Purpose: To generate isoparametric shape function data. (Single 

precision version of subroutine IHEXSD, section 4.27.8.43) 

3. Calling Sequence: CALL IHEXSS (ITYPE, SHP, DSHP, JINV, DETJ , EID, 

XI, ETA, ZETA, C0RD) 

See section 4.27.8.43 for description of arguments 
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4.46.8.39 Subroutine Name: SIHEX2 

1. Entry Point: SIHEX2 

2. Purpose: To perform final stress computations for the IHEX1, IHEX2, 

and IHEX elements. SIHEX2 must be called once for each 
stress point in an element. 

3. Calling Sequence: CALL SIHEX2 (TYPE, GPT) 

( 1 - IHEX1 ] 

TYPE = < 2 - IHEX2 > - integer-input 
( 3 - IHEX3 ) 

GPT - Array of grid point temperatures - real-input 

Section 4.46.10 Diagnostic Messages,, the following messages may additionally 
be issued by SDR2: 

2143, 4024 


89 



4.49 Functional Module DSMG1 (Differential Stiffness Generator - Phase 1) 

Section 4.49.7 Method, the first paragraph of this section should be replaced 

by the following: 

The module driver, DSMG1, is a very short routine whose only function 
is to call the two principal subroutines of the module, DS1 and DS1A, which 
accomplish the two phases of the module. The first phase of the module is 
incorporated in subroutine DS1. This routine creates the scratch file 
ECPTDS (GIN0 file number 301) by appending to each element in the ECPT 
data block for which differential stiffness is defined an element defor- 
mation, an average element loading temperature for all elements except the 
conical shell element and the isoparametric solid elements (a loading 
temperature at each grid point is appended for these elements) and the 
proper components of the displacement vector, UGV. It should be noted 
that although element deformations are defined only for rods, tubes, 
beams, and bars, an element deformation is attached to each element 
written on the ECPTDS scratch data block. The elements admissible to the 
ECPTDS scratch data block are: rods, beams, tubes, shear (but not twist) 

panels, triangular and quadrilateral elements (TRMEM's and QDMEM's), the 
combination membrane and plate triangular and quadrilateral elements 
(TRIA1, TRIA2, QUAD1, QUAD2) , the conical shell element, and the iso- 
parametric solid elements (IHEX1, IHEX2, and IHEX3) . 

Item 3 should read: 

3. If there is no temperature load, go to step 4. If there is a 
temperature load, the Grid Point Temperature Table data block, 

GPTT, is positioned to the proper thermal record. 

Sentence 10 of the first paragraph of item 6 should read: 

Then an element deformation number, an average element loading 
temperature (grid point temperatures for the conical shell and 
isoparametric solid elements) , and the displacement vector 
components are appended in core to the ECPT entry for the element. 
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The following two sections should be added: 
4.49.8.11 Subroutine Name: DIHEX 


1. Entry Point: DIHEX 

2. Purpose: To generate the element differential stiffness matrix 

for the isoparametric solid elements. 

3. Calling Sequence: CALL DIHEX (TYPE) 

( 1 - IHEX1 ) 

TYPE = < 2 - IHEX2 > - integer-input 
( 3 - IHEX3 ) 

4.49.8.12 Subroutine Name: IHEXSD 

See section 4.27.8.43 
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4.87 Structural Element Descriptions 


The following sections should be added: 

4.87.16 The IHEX1, IHEX2, and IHEX3 Elements 
4.87.16.1 Input Data for the IHEXi Elements 

1. The EST entries for IHEXi are: 


Symbol 


Description 

SIL^ , i=l,2, . 

. . ,NGP 

Scalar indices of the connected grid 
points. NGP = number of grid points 
per element. , 

N. 

l 


j 

X. 

" 1=1,2,. 


Referenced local coordinate system and 

. . ,NGP 

1 

location in basic coordinates of j 

Y. 

X 


connected grid points . 

Z. 

l 


1 

MID 


Material identification number 1 

CID 


Coordinate system identification number 
in which anisotropic material is defined. 

NIP 


i 

Number of integration points in each 1 

coordinate direction. J 

1 

AR 


1 

Maximum aspect ratio of element 

ALFA 


Maximum angle (in degrees) for face normals 

j 

BETA 


Maximum angle (in degrees) for mid-edge 
points 

GPT i , i=l,2, . 

. . ,NGP 

Grid point temperatures 
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2. Coordinate System Data 

The numbers N., X., Y., and Z. are used to calculate the three 
111 1 

by three global-to-basic coordinate transformation matrix [T.] 
th 

for the i — grid point via calls to either TRANSS or TRANSD. 

The CID is used to calculate the six by six global-to-basic 

M 

coordinate transformation matri [T ] for anisotropic material. 

If the CID references a cylindrical or spherical coordinate 
system, [T^T must be computed at each integration point. Other- 
wise, it need be computed once only for the element. 

3. Material Data 

Symbol Description 

[G ] Six by six stress-strain matrix defined 

m 

in coordinate system CID 
^ Mass density 

{a } Vector of six thermal expansion coefficients 

m 

defined in coordinate system CID 
T q Reference temperature 

g^ Structural damping coefficient 
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4.87.16.2 Basic Equations for IHEXi Elements 

1. Numerical Integration 

The method of Gaussian Quadrature is used to numerically integrate 

the element matrices. The weighting coefficients H_^ and abscissas 

s are listed in Table 1 at the conclusion of this Section 4.87.16 
i 

2. Element Coordinates 


The coordinates of the grid points in element coordinates are: 
For the IHEXI element: 


1 

2 

3 

4 

5 

6 

7 

8 


*i 

-1 

1 

1 

-1 

-1 

1 

1 

-1 


-1 

-1 

1 

1 

-1 

-1 

1 

1 


? i 

=1 

-1 

-1 

-1 

1 

1 

1 

1 
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For the IHEX2 element : 


i 


n i 


1 

-1 

-1 

-1 

2 

0 

-1 

-1 

3 

1 

-1 

-1 

4 

1 

0 

-1 

5 

1 

1 

-1 

6 

0 

1 

-1 

7 

-1 

1 

-1 

8 

-1 

0 

-1 

9 

-1 

-1 

0 

10 

1 

-1 

0 

11 

1 

1 

0 

12 

-1 

1 

0 

13 

-1 

-1 

1 

14 

0 

-1 

1 

15 

1 

-1 

1 

16 

1 

0 

1 

17 

1 

1 

1 

18 

0 

1 

1 

19 

-1 

1 

1 

20 

-1 

0 

1 
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For the IHEX3 element: 


i 


n i 


1' 

-1 

-1 

-1 

2 

-1/3 

-1 

-1 

3 

1/3 

-1 

-1 

4 

1 

-1 

-1 

5 

1 

-1/3 

-1 

6 

1 

1/3 

-1 

7 

1 

1 

-1 

8 

1/3 

1 

-1 

9 

-1/3 

1 

-1 

10 

-1 

1 

-1 

11 

-1 

1/3 

-1 

12 

-1 

-1/3 

-1 

13 

-1 

-1 

-1/3 

14 

1 

-1 

-1/3 

15 

1 

1 

-1/3 

16 

-1 

1 

-1/3 

17 

-1 

-1 

1/3 

18 

1 

; -i 

1/3 

19 

1 

i 

1/3 

20 

-1 

i 

1/3 

21 

-1 

-i 

1 

22 

-1/3 

-i 

1 

23 

1/3 

;-l 

1 

24 

1 

! -l 

1 

25 

1 

-1/3 

1 

26 

1 

1/3 

1 

27 

1 

"1 

1 

28 

1/3 

1 

1 

29 

-1/3 

1 

1 

30 

-1 

1 

1 

31 

-1 

1/3 

1 

32 

-1 

-1/3 

1 
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3. Interpolating Functions and Their Derivatives 

The interpolating function, or isoparametric shape function 

N. and the derivative of this function with respect to the 
i th 

element's coordinates is given for the i — grid point by the 
following tables: 2, 3, and 4. 

4. Jacobian Matrix 



The Jacobian matrix at any point (x,y,z) within the element is: 



Strain Displacement Relations 

th 

The transformations from displacements to strain for the i — 
grid point are: 
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where 



( 3 ) 


4.87.16.3 Stiffness and Mass Matrix Calculation for IHEXi Elements 

* 

The equation used in the stiffness matrix generation in global coordinates ‘is : 


[K. .] = [T.]' 

ij l 


NIP NIP NIP 

zee 

&=1 m=l n=l 


ii T M T M 

H. H H J [C.] L [T J 1 [G ] [T M ] [C . ] 
X. m n 1 1 l m 1 j 


. [T-] (4) 


The equation used in the coupled mass matrix generation in global coordinates is: 


[M..] = [T . ] 

ij 


NIP NIP NIP 

E E EpH, H m H n |J| 
36=1 m=l n=l 


N.N. * 0 0 

i 3 


0 N.N. 0 
i 3 


0 0 N.N. 

i J 


[T.] (5) 


] and [M „ ] are three by three partitions of the element matrices coupling 
the i*"* 1 and element grid points.* 


* No provision is made for the computation of a lumped mass matrix for the 
isoparametric hexahedron elements. Such a computation would necessitate 
an arbitrary distribution of mass reducing the accuracy advantages of 
these type elements. 
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4.87.16.4 Element Load Calculations for IHEXi Elements 


1. Thermal Force Vector (Subroutine IHEX of Module SSG1) 


The thermal force vector is generated with the following equation: 



NIP NIP NIP 


EEZ 

0=1 n =1 


H £ H m H n l J l t C i] T t TM ] T [G e ] fo e > * 




where {P.} is the load vector of length three for the i — element 
i th 

grid point and t ^ is the temperature at the j — element grid point. 

2. Pressure Force Vector (Subroutine PL0AD3 of Module SSGI) 

The generation of the pressure force vector is described in 
Section 4.41.8.9 

3. Heat Generation Vector for Heat Transfer Problems (Subroutine QIHEX 
of Module SSGI) 


The terms of the heat generation vector are given by: 

1 - 2 2 2 

: = -q LEE |j| n ( 7 ) 

SL=1 m=l n=l 

where Q is the external heat per unit volume. 
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4.87.16.5 Element Stress Calculations for IHEXi Elements 


1. Calculations performed in SIHEX1 (Phase 1 calculations) 

The Phase 1 calculations include the computation of stress matrices 
at each point within the element at which stresses are to be evaluated 

til 

The transformation matrix from displacements at the i — grid point to 
til 

stresses at the j — stress point is: 


[S j± ] = [TV [GJ [T M ] [C.] [T.] 


( 8 ) 


where [T . ] is the global to basic coordinate transformation matrix. 
All matrices on the right-hand side of this equation are computed at 
stress point j . 

The temperature to stress relation at the j— stress point is: 


{s > = -[TV [G ] {a } 

tj m m - 


(9) 


2. Calculations performed in SIHEX2 (Phase 2 calculations) 

The Phase 2 calculations include the computation of stresses at each 
stress point within the element. 
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w 


The^ equation for stress at the jrr; stress point is 


a 


NGP 


’ £ 'V ‘V + <V 


NGP 

£ N k (t k ' T o ) 
k=l k R ° 


(10) 


xy 


yz 


zx 


Where^U .} is the three by one global displacement vector at the i— 

S 1 th 

grid point and t^ is the loading temperature at the grid point. 

The three principal . stresses are the roots of the following cubic 
equation in S: 


S 3 - (a +a +a)S 2 +(ao + a a + aa - a 2 -a 2 -a ^ 2 ) S - 

x y z x y y z x z xy yz zx 


/ 


> 2 2 2 

Kaaa+ 2 aaa - a a - a a -aa ) = o 

x y z xy yz zx x yz y zx z xy 


The direction cosines of the normals to the principal plane on which 
each of the principal stresses is acting are found by solving: 


- (S fa x r o 

-a 

xy 

-a i 

zx 

-a 

xy 

(S - a y ) 

-a 

.yz 

-cr v 

*- zx 

-a 

yz 

(s -a z )_ 


using the constraint equation: 

2 2 2 
£ + m + n = 1 


J l 
m 


n 



( 12 ) 


( 13 ) 


( 11 ) 
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The mean stress or hydrostatic pressure is: 


a = - i (a +a + a )* 

n 3 x y z 


(14) 


The octahedral shear stress is: 


T o = + 0_) 2 + (S„ + C_) 2 + (S_ + a 


x n 


y n 




1/2 


(15) 


4.87.16.6' Differential Stiffness Calculations for IHEXi Elements 

The differential stiffness matrix in generalized coordinates is (the stresses 
are in the basic system) : 


[K u ] = 
w 


r m a 

-a 

-a 

yf j z 

xy 

xz 

-a 

a + a 

-a 

xy 

X Z 

o yz.. 

-a 

-a 

0+0 

xz 

yz 

X 


(16) 


The" transformation from displacements at the nodal points to rigid body 
rotation are: 


= -t 


0 


C 3i 


L -c 9 . 
£> 21 


" C 3i 


C li 


'2i 


- C li 


(17) 


where C, . , C„., and C„. are given in Section 4.87.16.2. 
li 2i 3i 

The three by three partitions of the global element differential stiffness 

th t h 

matrix couping the i — and j — element grid points are given by: 


[Ki/] = [T.] T 


2 2 2 

C'E I 

1=1 m=l n=l 


T, T, •£ iji [c, d ] T [k d ] [c d ] 


w 


[T j ] 


(18) 


* The minus sign is used to be consistent with the hydrostatic pressure 
computations for elements CTETRA, CWEDGE, CHEXA1, and CHEXA2. 
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4.87.16.7 Heat Transfer Calculations for IHEXi Elements 


The terms of the heat conductance matrix are given by: 


K. 


ij 


NIP NIP NIP 

= EEEh H„H„ |j| {C,} T [T M ] T [k] [T M ] {C'l 


£=1 m=l n=l 


m n 


( 19 ) 


where 


3N. 

i_ 

3x 


{C.} = 

1 


3N. 
• i 

3y 


( 20 ) 


3N. 

- i 




and [k] is the three by three matrix of thermal conductivity coefficients 
in global coordinates. ; 

The terms of the heat capacitance matrix are given by: 

_ J v 

Q- v O NIP NIP NIP 

B.. - Z jTj L, H.H H |j| H.cN. (21) 

l-l m-1 n=l 1 " n ^ J 


where c p is the thermal capacity per unit volume. 
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TABLE 1. GAUSSIAN QUADRATURE FORMULA 


111 n n n 

fff f(x, y, z)dx dy dz = ^ ^ X/ H f(s , s , s ) 

J -1 ±=1 k=l j=l J J K 1 

n 

Abscissa (s) 

Weight Coefficient (H) 

2 

±0.57735026919 

1.0 

3 

±0.77459666924 

0.0 

0.55555555555 

0.88888888888 

4 

±0.86113631159 

±0.33998104358 

0.34785484514 

0.65214515486 















TABLE 2. ISOPARAMETRIC SHAPE FUNCTIONS AND THEIR DERIVATIVES FOR THE IHEX1 ELEMENT - 8 GRID POINTS 



Where 



ISOPARAMETRIC SHAPE FUNCTIONS AND THEIR DERIVATIVES FOR THE IHEX2 ELEMENT - 20 GRID POINTS 



Where 







TABLE 4 . ISOPARAMETRIC SHAPE FUNCTIONS AND THEIR DERIVATIVES FOR THE IHEX3 ELEMENT - 32 GRID POINTS 
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INSTRUCTIONS FOR IMPLEMENTATION INTO NASTRAN LEVEL 16 


There are three major tasks involved in the implementation of the family 
of isoparametric solid elements into NASTRAN Level 16. These tasks are 
in addition to the normal implementation of subroutines and updates as 
delivered. Each of these tasks will involve changes in the source code. 
The first task is the installation of subroutine XIHEX for element matrix 
generation. The second is the installation of subroutine QIHEX for heat 
generation load computation. The third is the addition of the capability 
to generate element matrices and loads using general, 3-dimensional 
anisotropic materials, when such materials become available in NASTRAN. 

Subroutine XIHEX 


All program interfaces for this subroutine are based on MacNeal-Schwendler 
Corporation memo number DNH-4, dated June 21, 1972. The variables which 
control program branching to generate the various types of matrices are pre- 
sented in Table D.l. The GIN0 files which are passed to subroutine EMG0UT 
for writing the various matrices are shown in Table D.2. The work items 
which must be performed to implement XIHEX are listed below: 

1. Discard subroutines KIHEX , MIHEX, DIHEX, and TKTZTK . 

2. The local array GPTLD should contain the grid point temper- 
atures for computing stresses for the differential stiffness 
matrix. GPTLD must be initialized either by placing it in 
the argument list or in a common block. If GPTLD (1) = -1 
(integer) , the program assumes there is no temperature load. 

3. An integer pointer, UGV, must be initialized. UGV indicates 
the location of the first element of the single precision 
global displacement vector in open core. It is also used as 
the flag for computing the differential stiffness matrix (see 
Table D.l). UGV>0 to compute differential stiffness; UGV<0 
otherwise. 

4. The calling sequences for the four calls to EMG0UT should be 
checked. 

5. The matrices generated by XIHEX should be thoroughly checked 
against those generated by the Level 15 subroutines. 


Subroutine QIHEX 

No documentation on how this subroutine was to interface with the calling 
routine was available. Thus, the interfaces below were assumed. These 
must be modified as required to implement QIHEX. 

1. The element's EST entry is in common block ESTQ. 

2. The element load vector is generated and stored in array P. 

P is passed back to the calling routine in the argument list. 

The load is not inserted into a global load vector. 
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3. The load vector generated is for a unit value of Q, the external 
heat per unit volume. Therefore, the load vector passed to the 
calling routine must be multiplied by the actual value of Q. 


Anisotropic Materials 

When a 3-dimensional anisotropic material becomes available in NASTRAN, 
changes may be effected in the isoparametric solid element routines to make 
use of it. The following subroutines use the element material properties 
and thus would have to be changed: XIHEX, IHEX, and SIHEX1. In each of 

these subroutines, the logic for branching to computations for isotropic 
materials or anisotropic materials has already been implemented. However, 
certain changes and additions must still be made. First the interfaces with 
the material property utility PREMAT/MAT will have to be modified as re- 
quired. Second, new sections of code involving the actual computation with 
an anisotropic material must be added where indicated. The locations where 
new code must be added are indicated by comment cards with equal signs in 
columns 2-72. The computations which must’ be performed there are also 
indicated on comment cards. These computations are documented in the 
Theoretical Manual inserts in Appendix A and in the Programmer 1 s Manual 
inserts in Appendix C. The actual amount of code to be added is small and 
involves only the computation of small intermediate matrices using the 
anisotropic material properties. 
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TABLE D.l RULES FOR MATRIX GENERATION 


Rule 

Matrix Generated 

KGG>0 , UGV^O, HEAT=0 

Stiffness 

KGG>0 , UGV>0, HEAT=0 

Differential Stiffness 

KGG>0 , UGV^O, HEAT=1 

Conductance 

MGG>0 , UGV<0, HEAT=0 

Mass 

MGG>0, UGV<0 , HEAT=1 

Capacitance 

KGG>0, MGG>0 , UGV<0 , HEAT=0 

Stiffness and Mass 

KGG>0, MGG>0 , UGV<0, HEAT=1 

Conductance and Capacitance 


TABLE D. 2 MATRIX OUTPUT FILES 




'Matrix 

Output File 

Stiffness 

KGG 

Mass 

MGG 

Differential Stiffness 

KGG 

Conductance 

KGG 

Capacitance 

MGG 
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APPENDIX E 

DEMONSTRATION PROBLEM REPORTS 


Rigid Format - 1 Static Analysis 

Rigid Format - 3 Normal Modes 

Rigid Format - 5 Buckling Analysis 
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RIGID F0RMAT No. 1, Static Analysis 

Static Analysis of a Long Beam Using Hexahedron Isoparametric Elements 


A. Description 

This problem illustrates the use of NASTRAN's three types of hexahedron 
isoparametric solid elements (IHEX1, IHEX2, and IHEX3) to solve a long beam 
problem. The beam, shown in Figure 1, is fully constrained at the left end 
and free at the right end. Four separate loading conditions are considered 
in the analysis. In addition, solutions are calculated for each loading 
condition using three distinct finite element models of the beam. Each model 
is comprised of a different type of hexahedron isoparametric solid element 
(IHEX1, IHEX2 or IHEX3) . Results are presented for each of the three models, 
and compared with closed form solutions. 


Y 



Figure 1. Cantilever Beam 
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B . Input 


1. Parameters: 

L = 144 - Length 
D = 24 - Depth 

W = 12 - Width 

E = 30 x 10^ - Modulus of Elasticity 
y = .3 - Poisson's Ratio 

a = 1.428 x 10 - 5 - Coefficient of Thermal Expansion 
p = 7.535 x 10”^ - Mass Density 

2. Loads: 

a) Linear Thermal Gradient (Y-Direction) 

T = 120° @ Y = 0 , T = -120° @ Y = 24 

b) Uniform Temperature Rise 

AT = 100° 

c) Compressive Axial Pressure Load (Z-Direction) @ Z = 144 

P z = -42,837 

d) Transverse Pressure Load (Y-Direction) @ Y = 0 

P Y = 100 

3. Constraints: 

All grid points in the plane Z = 0 are constrained in the 
X, Y and Z directions. 

C. Modeling Techniques 

The finite element model for each type of isoparametric hexahedron 
element is shown in Figures 2, 3 and 4. A summary of the number of elements 
and grid points used in each model is given in Table 1 below. 


Table 1. Total Number of Elements and Grid Points for 
Each Type of Isoparametric Hexahedron Model 


Model Type 

Number of 
Elements 

Number of 
Grid Points 

Half-Band Width 
(Grid Points) 

Model 1 - IHEX1 Elements 

216 

364 

34 

Model 2 - IHEX2 Elements 

36 

275 

52 

Model 3 - IHEX3 Elements 

8 

148 

44 
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FIGURE 2. IHEX1 MODEL (216 ELEMENTS & 364 GRID POINTS) 
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FIGURE 3. IHEX2 MODEL (36 ELEMENTS & 275 GRID POINTS) 
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FIGURE 4. IHEX3 MODEL (8 ELEMENTS & 148 GRID POINTS) 
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D. Results 


Table 2 summarizes the maximum displacement obtained by NASTRAN for 
each of the three finite element models, and the approximate theoretical 
solution obtained by presented formulas. The results for the IHEX2 and 
IHEX3 models show a maximum error of 4 per cent of the theoretical solutions, 
and the IHEX1 model is within 8 per cent of the theoretical solutions. 


Table 2. Comparison of NASTRAN and Theoretical Solutions 


Loading Case 
Maximum Displacement 

NASTRAN 

- 

Model No. 1 
IHEXl 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

Theoretical* 

Solution 

Y-Disp. - Load Case a 

1.444 

1.548 

1.533 

1.481 

Z-Disp. - Load Case b 

.2113 

.2104 

.2088 

.2056 

Z-Disp. - Load Case c 

-.2039 

-.2042 

-.2047 

-.2056 

Y-Disp. - Load Case d 

. 1422 

.1561 

.1569 

.1586 


*Tlieoretical; Solutions 

Load Case a Load Case b Load Case c Load Case d 



6 y = 1.481 6 Z = .2056 6 Z = -.2056 6 y = .1586 
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RIGID FORMAT No . . 3 , Real. Eigenvalue Analysis 
Vibration Analysis of a. Long Beam. Using . Hexahedron Isoparametric Elements 

A. Description 

This problem illustrates . the use. of NASTRAN 's three. types of hexahedron 
isoparametric . solid . elements . (IHEX1, . IHEX2 and IHEX3) . to obtain the vibration 
modes . and frequencies . of . a cantilevered beam. The beam geometry and the three 
finite . element .‘models, employed, in the analysis are identical to the models 
used. in the previous . static analysis problem. 

B. Input . 

1. Parameters: 

L = 144 - Length 
D =24 - Depth 

W = 12 . — Width 

E = 30 x 10^ - Modulus of Elasticity 
y = .3 - Poisson's Ratio 

a = 1.428 x 10 - 5 _ Coefficient. of Thermal Expansion 
p = 7.535 x 10~4 - Mass Density 

2. Constraints: 

All grid points in the plane Z = 0 are constrained in the Z 
direction. In addition, the points on the line Z = 0, Y = 12 
are constrained in the Y direction, the points on the line Z = 0 
X = 0 are constrained in the X direction for the IHEX1 and IHEX3 
models, and the points on the line Z = 0, X = 6 are constrained 
in that direction for the IHEX2 model. 

3. Eigenvalue extraction data: 

Method: Inverse power 

Region of interest: 0 < f < 120 

Number of desired roots: 3 

C. Results 

Table 1 and Figures 1 and 2 summarize the NASTRAN results for the three 
lowest natural frequencies and associated mode shapes. Table 1 lists the 
theoretical frequencies and the NASTRAN frequencies for each type of hexa- 
hedron isoparametric element model. The results for IHEX2 and IHEX3 elements 
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are within 3 per cent of the theoretical solution. The results for the 
IHEX1 elements are within 4 per cent of the theoretical solution for bend- 
ing in the Y direction, but overestimate the first and second mode frequen- 
cies for bending in the X direction by 14 and 17 per cent, respectively. 

The latter is probably due to not having sufficient grid points through the 
beam in the X direction. This problem demonstrates that the IHEX2 and IHEX3 
elements are better choices for beam bending problems than the IHEX1 element. 

Plots of the first three beam bending modes are shown in Figures 1 and 
2 for each type of finite element model, and are compared with the theoret- 
ical solution. All NASTRAN results show excellent correlation. 


Table 1. Comparison of Natural Frequencies from NASTRAN 

with Theoretical 






NASTRAN 




Mode 

No. 

Description 

Model No. 1 
IHEX1 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

Theoretical 

Solution* 

(cps) 



Freq. 

% 

Error 

Freq. 

% 

Error 

Freq. 

% 

Error 

1 

First Bending Mode 
in the X-Direction 

22.0 

+18, 

18.6 

0 

18.6 

0 

18.6 

2 

First Bending Mode 
in the Y-Direction 

38.3 

+4 

36.5 

-2 

36.5 

-2 

37.3 

3 

Second Bending 
Mode in the 
X-Direction 

135.3 

+16 

114.3 

-2 

113.3 

-3 

116.8 


~W. C. Hurty and M. F. Rubenstein, Dynamics of Structures , Prentice-Hall, 
Inc., Englewood Cliffs, New Jersey, 1964, (Chapter 5). 
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FIGURE 1. FIRST BEAM VIBRATION BENDING MODES, X AND Y DIRECTIONS 



IHEX1 Model 
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FIGURE 2. SECOND BEAM VIBRATION BENDING MODE, X DIRECTION 


RIGID FORMAT No. 5, Buckling Analysis 


Buckling Analysis of a. Long Beam. Using. Hexahedron Isoparametric Elements 


A. Description 

This problem illustrates the use. of NASTRAN's three types of hexahedron 
isoparametric. solid elements (IHEX1, IHEX2 and IHEX3) to obtain the axial 
buckling load of a. cantilevered beam. The beam geometry and the three finite 
element models employed. in the analysis. are identical to the models used in the 
previous static . analysis problem. 


B. - Input '> - - ... - 

1. Parameters: *• - _ 

L = 144 - Length 

D = 24 - Depth 
W = 12 - Width 

E - 30 x 10^ - Modulus. of Elasticity 
y.= .3 - Poisson's Ratio 

a.- 1.428 x 10 -5 - Coefficient of Thermal Expansion 
p = 7.535 x 10 - ^ - Mass Density 

2 . iLoads': . '■ 

Compressive Axial Pressure Load. (Z-Direction) @ Z = 144 
P z = 42,837 psi 

3. Constraints: 

All grid points in the plane Z = 0 are constrained in the Z direction. 
In addition the points on. the line Z = 0, Y = 12 are constrained in the 
Y - direction, the points on the line Z = 0, X = 0 are constrained in 
the X - direction for the IHEXl and IHEX3 models and the points on 
the line Z = 0, X = 6 are constrained in that direction for the IHEX2 
model. 

4. Eigenvalue extraction data: 

Method: Inverse Power 
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Region of Interest: 0 < A < 10 

Number of desired roots: 3 

C . Results 

The critical buckling load* for a long cantilevered beam is given by 

„ , , . , , m2 tt' 2 El /, \ 

Buckling Load = — ^ (1) 

Substituting the appropriate . parameters into (1) yield: 

For m = 1, Buckling ..Load (X - Dir.) = 1.234 x 10^ lbs. 

For m = 1, Buckling Load (Y - Dir^) = 4.955 x 10 A lbs. 

For m = 2, Buckling Load (X - Dir.) = 1.110 x 10^ lbs. 

y ■ or First Mode P cr i t i ca i ( x ” Dir *) = 42837 psi 

First Mode. ? z critical ( Y " Dir.) = 4(42837) = 170,948 psi 

Second Mode P , (X - Dir.) = 9(42837) = 385,533 psi 

z critical 

Since the compressive pressure . load applied to the beam is equal to the 
minimum buckling load, the eigenvalue . which NASTRAN should obtain is 
A = 1.0 for buckling in the X-direction, A = 4.0 for buckling in the 
Y-direction, and A = 9.0 for buckling in the second mode, X-direction. 

Table 1 and Figure 1 and 2 summarize the NASTRAN results for buckling in 
the X and.Y directions. Table 1 lists the NASTRAN calculated buckling parameter, 
A for each type of isoparametric element model. The IHEX2 and IHEX3 element 
results are excellent and are within .7 percent of the theoretical solution. 
However, as in the vibration problem, the IHEX 1 element model overestimates the 
buckling load in the Y direction by 10 percent, and the buckling load in the 
X-direction by 41 and 42 percent respectively. Again, the latter is most 
likely due to the use of insufficient grid points through the X-direction of the 
beam. 


* S.P. Timoshenko and J. M. Gere, Theory of Elastic Stability , McGraw-Hill 
Book Company, Inc., New York, 1961 (Chapter 2). 
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The fundamental buckling modes in the X and Y directions are compared 
with the theoretical in Figure 1. The second buckling mode in the X direction 
is shown in Figure 2. The correlation with theory is excellent. 


Table 1. Comparison of Buckling Loads from NASTRAN 
with Theoretical 


. 

Mode 

No. 

Description 

NASTRAN 

Theoretical 

Solution 

A 

Model No. 1 
IHEX1 
Elements 

Model No. 2 
IHEX2 
Elements 

Model No. 3 
IHEX3 
Elements 

m 

% 

Error 

A 

% 

Error 

A 

% 

Error 

■ 

. 

X-Direction 

1.406 

40.6 

1.002 

, .2 

1.001 

.1 

1.0 

2 

Y-Direction 

4.391 

kO 
; 00 

3.981 

.5 

3.979 

.5 

4.0 

3 

X-Direction 

12 . 809 

42.3 

9.037 

.4 

8.934 

.7 

9.0 
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FIGURE 1. FIRST BEAM BUCKLING BENDING MODE, X AND Y DIRECTIONS 




129 


FIGURE 2. SECOND BEAM BUCKLING BENDING MODE, X DIRECTION 



